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Summary 


This  document  presents  the  work  performed  under  the  auspices  of  the  Air  Force  Office  of  Scientific 
Research  (AFOSR)  through  the  European  Office  of  Aerospace  Research  and  Development  (EOARD)  with 
award  number  FA8655-10-1-3005.  The  project  focused  on  the  mixing  behavior  of  fluids  under  transcritical 
flow  conditions.  This  topic  is  relevant  to  rocket  propulsion  applications  and  it  is  also  applicable  to 
other  high  pressure  systems  such  as  gas  turbine  combustors.  Extensive  work  in  this  general  area  of 
propulsion  has  been  performed  at  the  EM2C  (Energetique,  Moleculaire,  Macroscopique  et  Combustion) 
laboratory  at  Ecole  Centrale  Paris  (ECP)  which  has  suitable  complemented  the  experimental  work  that 
has  been  recently  conducted  at  the  Air  Force  Research  Lab  (AFRL)  at  Edwards,  Air  Force  Base  (AFB) 
in  California.  The  report,  which  is  divided  in  three  major  parts,  encompasses  the  results  of  modeling 
efforts,  further  analysis  of  experimental  data,  and  comparisons  between  the  two.  The  first  part  deals  with 
the  Large-Eddy  Simulation  (LES)  of  non-reactive  shear  coaxial  jets,  the  second  part  expands  on  the  area 
of  spectral  analysis  and  the  final  part  presents  the  work  performed  using  an  in-house  code  developed  by 
the  Numerical  Combustion  Group  of  the  EM2C  laboratory.  This  code  is  based  on  the  Direct  Numerical 
Simulation  (DNS)  technique. 
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Introduction 


The  following  report  is  divided  into  three  parts  which  present  different  aspects  of  non-reactive  tran- 
scritical  flows  with  applications  in  liquid  rocket  engine  environments  and  other  propulsion  and  power 
generation  systems.  Part  I  of  the  report  is  concerned  with  the  large  eddy  simulation  of  coaxial  jets  in¬ 
jected  under  transcritical  conditions.  In  particular,  three  transcritical  experimental  cases  performed  by 
the  Combustion  Devices  Group  at  AFRL  were  modeled.  This  type  of  flows  feature  a  shear  coaxial  geom¬ 
etry  that  presents  many  interesting  challenges  such  as  understanding  and  quantifying  mixing  behavior 
and  spreading  of  the  inner  jet  into  the  outer  jet  and  the  surrounding  fluid.  The  thermodynamics  of  flows 
near  the  critical  point  of  the  fluid  -  a  point  which  indicates  the  beginning  of  a  thermodynamic  region 
where  the  state  of  the  fluid  cannot  be  distinguished  between  liquid  and  gas  -  are  also  complex,  with  some 
properties  acquiring  unusual  values  such  as  zero  surface  tension  and  zero  heat  of  vaporization.  The  flow 
dynamics  of  such  cryogenic,  transcritical  shear  coaxial  jets  are  of  interest  because  of  their  widespread 
use  as  propellant  injectors.  The  latter  also  makes  them  ideal  candidates  for  the  analysis  of  combustion 
instabilities.  It  is  desired  to  have  a  better  understanding  of  how  the  inner  and  outer  streams  of  a  coaxial 
jet  interact,  mix  and  burn  in  such  conditions  to  be  able  to  model  spatial  and  temporal  release  of  energy 
more  accurately.  This  in  turn  will  aid  in  analyzing  coupling  phenomena  between  the  injector  flow  and  the 
combustion  chamber  modes  which  are  mechanisms  that  produce  oscillations  in  pressure  and  temperature 
that  could  lead  to  undesired  phenomena  such  as  very  high  heat  flux  at  the  combustion  chamber  walls. 
Numerical  simulation  helps  uncover  details  that  cannot  be  found  using  experimental  techniques,  such  as 
knowing  the  value  of  the  velocity  vector  at  any  desired  position  in  the  flow  and  easily  extracting  spectral 
properties  from  a  given  location. 

These  spectral  properties  are  analyzed  in  Parts  I  and  II  of  this  report.  The  spectral  analysis  of  Part 
I  focuses  on  the  analysis  of  the  flow  structures  and  their  associated  spectral  content  as  observed  in  the 
simulations.  For  cases  with  no  acoustic  modulation,  numerical  probes  were  placed  in  points  of  interest 
within  the  flow  to  analyze  their  spectral  content.  For  instance,  in  regions  where  a  strong  recirculation 
zone  is  observed  or  at  a  point  where  a  mixing  layer  is  present.  Using  this  analysis,  preferred  modes 
or  frequencies  of  the  system  are  retrieved  and  associated  Strouhal  numbers  calculated.  The  spectral 
analysis  in  Part  I  also  includes  simulations  with  acoustic  modulation  for  the  purpose  of  understanding 
more  the  combustion  instability  phenomena  as  described  above.  For  the  cases  with  acoustic  modulation, 
the  spectral  content  was  obtained  not  only  for  one  point  but  also  for  an  entire  line  of  locations,  in  order 
to  follow  the  dynamics  of  the  flow  along  that  particular  axis.  Results  are  shown  at  the  end  of  Part  I  of 
this  report.  All  contents  of  Part  I  were  submitted  as  an  article  for  publication  in  the  journal  Physics  of 
Fluids  and  are  currently  under  review. 

Part  II  of  this  report  deals  with  the  analysis  of  the  spectral  content  of  coaxial  jets  injected  under 
transcritical  conditions  with  the  same  goal  of  identifying  any  periodicities  in  our  data  sets.  The  data  sets 
analyzed  now  include  some  of  the  experimental  data  gathered  at  AFRL.  The  signal  from  which  spectral 
content  is  extracted  is  the  intensity  of  light  at  each  pixel  from  images  captured  by  a  high-speed  camera. 
Spectral  maps  are  obtained  from  experimental  cases  with  and  without  acoustic  modulation.  The  same 
axes  as  those  used  in  Part  I  of  this  report  were  used  to  allow  comparison.  Spectral  analysis  is  a  technique 
that  could  be  used  extensively  with  other  experimental  datasets  already  obtained  at  AFRL  to  analyze 
frequency  content.  The  results  of  this  technique  applied  to  cases  without  acoustic  modulation  and  their 
comparison  with  selected  simulation  results  are  shown  first.  The  second  set  of  results  consists  of  cases 
with  acoustic  modulation.  Spectral  maps  of  both  experimental  data  and  simulation  results  are  shown. 
It  is  found  that  both  techniques  clearly  identify  the  modulation  frequency.  As  an  interesting  note,  given 
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the  difference  in  the  nature  of  the  signal  used  from  each  type  of  data,  the  quantitative  comparison  of  the 
physical  location  of  the  mode  for  one  of  the  cases  is  remarkable. 

The  last  part  of  this  report  focuses  on  the  use  of  a  direct  numerical  simulation  tool  to  model  nitrogen 
flows.  DNS  is  often  used  to  evaluate  certain  flow  properties  from  which  suitable  subgrid-scale  (SGS) 
models  are  developed  to  be  used  in  LES.  These  “subgrid-scale”  models  are  needed  in  LES  because  an 
LES  grid  is  too  coarse  to  resolve  the  smallest  eddies  and  flow  patterns  present  in  the  flow.  Since  the 
way  in  which  the  flow  behaves  and  turns  at  the  smallest  scales  is  responsible  for  an  important  energy 
transfer  step  (of  kinetic  energy  into  heat)  within  the  flow  the  need  to  know  how  these  smallest  eddies 
and  flow  patterns  behave  is  very  important.  DNS  grids  do  resolve  even  the  smallest  scales  of  the  flow 
so  they  are  able  to  take  into  account  this  energy  transfer  within  the  flow;  however,  using  a  DNS  grid  to 
model  the  flow  of  a  rocket  engine  combustion  chamber  or  around  an  aircraft  engine  turbine  blade  would 
be  impossible  given  the  number  of  grid  points  such  a  problem  would  require.  This  is  why  in  practice 
DNS  is  limited  to  very  simple  geometrical  domains. 

Nonetheless,  LES  of  transcritical  and  supercritical  flows  is  feasible  and  DNS  has  been  used  to  obtain 
SGS  models  to  incorporate  them  to  LES  codes.  Parts  I  and  III  discuss  some  of  the  groups  involved  in 
this  type  of  research.  In  Part  III,  the  direct  simulation  of  nitrogen  flow  using  a  DNS  tool  developed  at  the 
EM2C  laboratory  of  ECP  is  considered.  It  explains  its  capabilities,  the  different  thermodynamic  modules 
available  for  the  user,  and  the  validation  performed  using  a  one-dimensional  direct  numerical  simulation. 
This  part  contains  results  of  a  2D  simulation  of  the  evolution  of  a  shear  layer  formed  between  two  nitrogen 
streams  flowing  in  opposite  directions.  The  simulation  provides  detailed  evolution  of  a  smooth  transition 
between  the  two  flows  into  large  vortices  which  eventually  merge  to  form  larger  ones  highlighting  the 
mechanics  of  vortex  pairing. 
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Part  I 

Modeling  of  Transcritical  Mixing 
using  Large  Eddy  Simulation 
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Experiments  and  numerical  simulation  of  mixing 
under  supercritical  conditions 

Authors:  T.  Schmitt1,  J.  Rodriguez1,2,  I.  A.  Leyva2  and  S.  Candel1 

XEM2C,  CNRS,  Ecole  Centrale  Paris,  92295  Chatenay-Malabry,  France, 
2AFRL/RZSA,  Edwards  AFB,  CA.  93524,  USA 

Article  submitted  to  the  journal  Physics  of  Fluids. 
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Chapter  1 

Abstract 


Supercritical  conditions  designate  a  situation  where  the  working  fluid  pressure  is  above  the  critical  point. 
Among  these  conditions,  it  is  interesting  to  identify  a  transcritical  range  which  corresponds  to  cases  where 
the  pressure  is  above  the  critical  point,  but  the  injection  temperature  is  below  the  critical  value.  This 
situation  is  of  special  interest  because  it  raises  fundamental  issues  which  have  technological  relevance 
in  the  analysis  of  flows  in  liquid  rocket  engines.  This  situation  is  here  envisaged  by  analyzing  the 
behavior  of  a  nitrogen  shear  coaxial  jet  comprising  an  inner  stream  injected  at  temperatures  close  to 
the  critical  temperature  and  a  coaxial  flow  at  a  higher  temperature.  Experiments  are  carried  out  both 
in  the  absence  of  external  modulation  and  by  imposing  a  large  amplitude  transverse  acoustic  field.  Real 
gas  Large  Eddy  Simulations  are  performed  for  selected  experiments.  The  combination  of  experiments 
and  calculations  is  used  to  evaluate  effects  of  injector  geometry  and  operating  parameters.  Calculations 
retrieve  what  is  observed  experimentally  when  the  momentum  flux  ratio  of  the  outer  to  the  inner  stream 
J  =  (peul) / (piuf)  is  varied.  Results  exhibit  the  change  in  flow  structure  and  the  development  of  a 
recirculation  region  when  this  parameter  exceeds  a  critical  value.  The  instantaneous  flow  patterns  for 
different  momentum  flux  ratios  are  used  in  a  second  stage  to  characterize  the  dynamical  behavior  of 
the  flow  in  terms  of  power  spectral  density  of  velocity  and  density  fluctuations.  Results  obtained  under 
acoustic  modulation  provide  insight  on  mixing  enhancement  of  coaxial  streams  with  a  view  of  its  possible 
consequences  in  high  frequency  combustion  instabilities.  It  is  shown  in  particular  that  the  presence  of 
strong  acoustic  modulations  notably  reduces  the  high  density  jet  core  length,  indicating  an  increased 
mixing  efficiency.  This  behavior  is  more  pronounced  when  the  jet  is  placed  at  the  location  of  maximum 
transverse  velocity  fluctuations. 
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Chapter  2 

Introduction 


Flows  formed  by  coaxial  injectors  raise  many  fundamental  issues.  Such  flows  are  also  of  practical  relevance 
for  high  performance  cryogenic  liquid  rocket  engines  (LREs).  It  is  important  to  understand  the  basic 
processes  taking  place  in  these  flows  to  improve  the  design  methodology  for  future  cryogenic  LREs. 
Another  potential  application  of  research  reported  in  this  article  lies  in  the  area  of  combustion  instability. 
It  is  known  that  mechanisms  driving  combustion  oscillations  result  from  the  coupling  between  the  injector 
flow,  the  combustion  process  and  the  chamber  acoustic  modes.  The  manner  in  which  the  inner  and 
outer  streams  originating  from  a  coaxial  element  interact  with  each  other  and  with  their  surroundings 
may  affect  the  temporal  and  spatial  release  of  chemical  energy  within  the  chamber  leading  to  undesired 
oscillations  associated  with  a  resonant  acoustic  motion.  In  most  LREs,  performance  has  been  enhanced 
by  augmenting  the  chamber  pressure.  Current  LOx/H2  engines  like  the  SSME,  Vulcain2  and  RS-68 
operate  at  values  exceeding  the  critical  pressure  of  their  individual  propellants.  The  inner  jet,  oxygen,  is 
generally  injected  in  a  transcritical  or  liquid-like  state  (z.e.  at  a  pressure  higher  than  the  critical  value  but 
at  a  temperature  lower  than  the  critical  temperature)  and  it  is  surrounded  by  a  higher- velocity  hydrogen 
stream  at  supercritical  temperature  [72].  The  inner  fluid  stream  injected  in  a  transcritical  form  evolves 
to  a  supercritical  state  as  its  temperature  rises  because  of  mixing  and  possibly  combustion.  Under  such 
thermodynamic  conditions,  the  system  exhibits  particular  features  that  differ  from  those  of  a  two-phase 
flow  or  of  a  pure  gaseous  injection.  The  objective  of  the  present  work  is  to  compare  experiments  and 
simulations  in  which  shear  coaxial  jet  configurations  are  established  at  pressures  above  the  critical  point 
of  the  injected  fluid,  but  at  transcritical  temperatures  without  and  with  acoustic  excitation.  For  this 
study,  nitrogen  (N2)  is  the  only  fluid  used  in  the  inner,  outer  and  chamber  flows. 

At  this  point  it  is  worth  reviewing  previous  work  dealing  with  supercritical  injection.  Starting  with 
single  round-jets,  experiments  have  been  carried  out  at  AFRL  (Air  Force  Research  Laboratory)  [9]  and 
DLR  (Deutsches  Zentrum  fiir  Luft-  und  Raumfahrt)  [29,  41]  (see  also  Oschwald  et  al  [43]  for  a  review  up 
to  2006).  Some  recent  investigations  are  reported  by  Segal  and  Polikov  [61]  and  Roy  and  Segal  [56].  These 
studies  indicate  that  the  liquid  break-up  and  atomization  mechanism  which  prevail  at  subcritical  pressure 
is  no  longer  observable  because  surface  tension  and  latent  heat  of  vaporization  vanish  at  supercritical 
pressures  [47].  The  jet  of  fluid  then  dissolves  in  the  ambient  gas,  with  no  evidence  of  droplets.  The  flow 
features  “comb-like”  structures  at  the  edge  of  the  dense  stream  a  type  of  pattern  which  is  not  observed  at 
subcritical  pressures  [9].  The  geometry  of  the  flow  is  reminiscent  of  that  of  a  variable  density  turbulent 
gas  stream  [9].  This  behavior  has  been  confirmed  by  experimental  measurements  of  visual  spreading 
rates  in  the  initial  region  which  are  found  to  be  consistent  with  theoretical  growth  rate  expressions  for 
incompressible  variable-density  turbulent  mixing  layers.  Quantitative  measurements  of  density  [4,  41,  29], 
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spreading  rate  and  axial  decay  coefficients  obtained  by  DLR  [43,  40]  using  spontaneous  Raman  scattering 
also  convey  the  same  type  of  information. 

Coaxial  injection  is  investigated  by  several  researchers  at  AFRL  [14,  13,  25,  24,  55,  54]  who  exam¬ 
ine  the  single-species  case,  and  by  Mayer  and  Smith  [30]  who  examine  a  stream  of  nitrogen  injected 
at  low  temperature  (high  density)  surrounded  by  a  high  velocity  stream  of  moderate  temperature  he¬ 
lium.  Quantitative  measurement  of  species  density  are  reported  by  Oschwald  et  al  [42]  for  a  coaxial 
nitrogen  /  hydrogen  injection,  using  spontaneous  Raman  Scattering.  When  compared  with  the  round 
jet  configuration,  a  strong  reduction  of  the  jet  density  potential  core  length  is  observed  [43].  For  shear 
coaxial  jets,  it  has  been  found  that  the  momentum  flux  ratio  between  the  outer  stream  and  the  inner  jet 
(J  =  (peul) / (piuf))  is  of  major  importance  for  the  mixing  efficiency  [22,  69].  It  has  also  been  found  that 
combustion  is  more  stable  and  more  efficient  at  high  values  of  this  quantity  [5].  Different  flow  regimes 
may  be  obtained  by  varying  this  ratio  [51]  and  that  the  inner  jet  core  length  decreases  as  the  momentum 
flux  ratio  increases  [13,  25]. 

The  dynamics  of  coaxial  transcritical  jets  is  of  interest  for  the  study  of  combustion  instability.  This  can 
be  examined  by  considering  the  interaction  of  acoustic  waves  with  transcritical  flows  and  its  influence  on 
mixing  efficiency.  Effects  of  an  acoustic  modulation  on  a  single  round  jet  at  sub-  and  supercritical  pressure 
is  investigated  by  Chehroudi  and  Talley  [8]  who  find  that  effects  on  the  jet  structure  are  more  profound 
in  the  subcritical  case.  The  impact  of  the  modulation  is  reduced  as  the  initial  jet  velocity  is  augmented. 
The  impact  of  acoustics  on  a  non-reacting  shear  coaxial  jet  has  been  investigated  experimentally  by  Leyva 
and  Rodriguez  [25,  24,  55,  54].  They  have  assessed  the  effects  of  acoustics  on  the  jet  by  measuring  the 
spreading  angle  of  the  inner  and  outer  jets  and  the  length  of  the  inner  jet  before  its  first  break-up.  It  has 
been  found  that  the  injector  geometry  has  an  effect  on  the  susceptibility  of  the  coaxial  flows  to  acoustic 
modulation.  For  the  injector  geometry  studied  here,  acoustics  have  a  measurable  effect  on  the  reduction 
of  the  inner  jet  dark  core  length  for  an  intermediate  value  of  J  at  pressures  varying  from  1.5  to  5.0  MPa. 
Reactive  cases  were  also  experimentally  investigated  (see  [5,  19,  67,  64,  65]  for  example,  and  [53,  32]  for 
acoustically  modulated  cases)  but  will  not  be  detailed  here  since  they  are  out  of  the  scope  of  the  present 
article. 

The  numerical  modeling  of  such  complex  flows  has  been  considered  more  recently  by  various  groups. 
It  is  first  noted  that  in  the  transcritical  range,  thermodynamic  properties  notably  differ  from  those  of  a 
perfect  gas  and  cannot  be  accurately  represented  with  the  standard  perfect  gas  (PG)  equation  of  state 
(EOS).  Thermodynamics  and  “real-gas”  (RG)  EOS  are  required  [47,  15,  71].  Cubic  equations  of  state 
are  generally  used  to  this  purpose  [37],  but  more  precise  (and  numerically  expensive)  equations  may  also 
be  employed  [71].  This  issue  and  the  modeling  of  transport  properties  are  extensively  reviewed  in  [1], 
[37],  [31]  and  [2]. 

The  models  described  in  the  previous  studies  are  used  in  various  calculations.  Mixing  at  supercritical 
pressure  is  studied  with  direct  numerical  simulation  in  a  mixing  layer  configuration,  to  show  the  role  of 
density  gradients  on  the  global  layer  stability  and  turbulence  characteristics  and  the  strong  effect  on  the 
molecular  mixing  of  the  reduced  species  diffusion  near  critical  conditions  [2] .  This  highlights  the  need  for 
approriate  subgrid-scale  models  for  Large  Eddy  Simulation  (LES)  under  supercritical  pressure  conditions 
[63,  68]. 

A  few  large  eddy  simulations  (LES)  have  been  performed  under  transcritical  conditions.  An  extensive 
review  of  non-reacting  LES  in  the  period  ending  in  2006  can  be  found  in  Zong  et  al.  [76].  A  single  nitrogen 
round  jet  was  studied  by  Zong  [74,  73]  and  recently  by  Schmitt  et  al  [59].  Again,  the  stabilizing  effect  of 
the  density  gradient,  and  its  role  on  turbulent  energy  redistribution  along  the  mixing  layer  was  identified 
[74].  Shear  coaxial  injection  of  oxygen  and  methane  has  also  been  simulated  by  Zong  and  Yang  [75].  This 
has  been  continued  by  considering  a  non-reactive  coaxial  injector  submitted  to  acoustic  modulations  [26] 
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where  a  significant  reduction  of  the  jet  intact  core  length  is  observed  when  acoustic  perturbations  are 
present. 

One  also  notes  that  LES  have  been  used  for  reactive  cases  at  supercritical  pressures.  Shear  coaxial  jet 
flames  were  investigated  by  Oefelein  [37,  38],  Oefelein  and  Yang  [36],  Zong  and  Yang  [77],  Matsuyama 
et  al  [28],  Masquelet  et  al  [27]  and  Schmitt  et  al  [57].  The  problem  is  also  considered  in  the  Reynolds 
Average  Navier- Stokes  (RANS)  framework  (see  for  example  [16],  [11]  or  [49]). 

The  present  investigation  focuses  on  the  mixing  behavior  in  the  flow  originating  from  a  shear  coaxial 
injector  in  the  absence  or  in  the  presence  of  external  acoustic  modulation.  Three  operating  points  are 
considered  corresponding  to  different  momentum  flux  ratios.  The  coaxial  jet  geometry  corresponds  to  a 
thick  inner  jet  post  (the  wall  which  separates  the  inner  and  outer  streams  in  the  coaxial  configuration  has 
a  thickness  of  the  order  of  the  jet  diameter).  In  this  article  each  case  is  defined  by  the  chamber  pressure, 
the  momentum  flux  ratio,  and  the  coaxial  jet  geometry  which  is  fixed  in  the  present  study.  The  main 
components  of  the  experimental  facility  are  reviewed  in  chapter  3.  Balance  equations  and  models  used  in 
the  simulations  are  briefly  described  in  chapter  4.  The  definition  of  the  computational  domain  and  other 
numerical  aspects  are  considered  in  chapter  5.  The  influence  of  the  momentum  flux  ratio  is  examined 
in  chapter  6.  Results  obtained  are  used  to  examine  the  dynamics  of  the  flow  and  identify  the  spectral 
content  in  the  various  regions.  Effects  of  a  transverse  acoustic  modulation  are  investigated  in  chapter  7. 
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Chapter  3 


Experimental  configuration 


The  experimental  work  described  in  this  article  was  carried  out  at  the  Supercritical  Cryogenic  Test 
Facility  (EC-4)  of  the  Air  Force  Research  Laboratory  at  Edwards  Air  Force  Base,  California.  The  facility 
features  a  test  chamber  that  can  achieve  pressures  over  5  MPa  and  cryogenic  injection  at  temperatures 
below  120  K.  The  photograph  in  Fig.  3.1(a)  shows  the  main  test  chamber  and  several  key  components 
of  the  facility.  A  schematic  view  of  the  test  chamber  is  provided  in  Fig.  3.1(b).  The  facility  was  used 
to  obtain  a  shear  coaxial  jet  flow  where  the  inner  jet  was  cooled  down  to  be  denser  than  the  outer  jet 
under  typical  combustion  chamber  pressures,  resembling  conditions  found  in  practice  where  the  inner 
jet  oxidizer  is  injected  at  a  significantly  lower  temperature  compared  to  that  of  the  outer  jet  stream.  A 
second  equally  important  objective  was  to  include  the  capability  to  perform  acoustic  forcing.  When  this 
was  achieved  the  coaxial  jet  flow  parameters  were  varied  and  the  coaxial  jet  behavior  was  observed  with 
and  without  acoustic  excitation. 


Chamber  vent 


Figure  3.1:  (a)  Photograph  of  the  main  test  chamber  and  key  instrumentation  of  EC4,  the  supercritical 
cryogenic  facility  at  AFRL,  Edwards  AFB.  (b)  Flow  diagram  of  the  test  chamber  and  its  vicinity. 


To  run  a  test,  a  laboratory- wide  high  pressure  gaseous  nitrogen  supply  was  fed  into  the  inner  and 
outer  flow  lines  of  the  coaxial  jet  and  the  flow  rates  adjusted  to  correspond  to  those  of  the  desired  test 
conditions.  The  mass  flow  rates  through  the  inner  and  outer  jets  were  measured  with  Porter  mass  flow 
meters  (models  122  and  123-DKASVDAA).  The  same  gaseous  nitrogen  supply  used  for  the  coaxial  jet 
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was  also  employed  to  increase  the  pressure  in  the  test  chamber  and  achieve  the  required  level,  which  was 
measured  with  a  Stellar  1500  transducer.  The  jet  temperatures  were  controlled  by  adjusting  the  flow 
rates  of  liquid  nitrogen  through  the  heat  exchangers  placed  on  the  feed  lines.  Densities  of  the  two  jets 
were  deduced  from  the  measured  temperatures  and  pressure  using  the  National  Institute  of  Standards  and 
Technology’s  thermophysical  properties  of  fluid  systems  online  database  [23].  These  values  in  conjunction 
with  the  mass  flow  rate  data  provided  the  injection  velocities  and  the  corresponding  Reynolds  numbers, 
velocity  and  momentum  flux  ratios.  Several  thermocouples  were  used  across  the  heat  exchanger  and  in 
other  locations  to  keep  track  of  the  conditions  of  the  flow  in  order  to  maintain  the  required  flow  properties 
for  each  test  case.  Appendix  A  provides  a  summary  of  the  operating  conditions  achieved  for  all  cases 
reported  in  the  present  article.  Two  images  of  the  shear  coaxial  injector  can  be  seen  in  Fig.  3.2.  The 
first,  Fig.  3.2(a),  gives  a  full  view  before  assembly  and  the  second,  Fig.  3.2(b)  shows  a  plane  view  of  the 
end  section.  The  inner  jet  diameter  is  d{  =  0.51  mm.  The  outer  jet  has  an  inner  diameter,  of  1.59  mm 
and  outer  diameter  of  2.42  mm.  The  injector  outer  diameter  was  3.18  mm.  The  length  to  inner  diameter 
ratio  was  100  for  the  inner  jet  and  34  for  the  outer  jet  based  on  the  hydraulic  diameter.  The  coaxial 
injector  was  installed  so  that  the  inner  and  outer  jets  were  nominally  concentric.  The  inner  jet  exit  plane 
was  recessed  by  0.3  mm  from  the  outer  jet  exit  plane.  This  recess  length  was  chosen  to  mimic  realistic 
coaxial  jet  configurations  used  in  practice. 


Figure  3.2:  (a)  Photograph  of  the  shear  coaxial  injector  used  in  this  study,  (b)  Exit  plane  view  of  the 
injector. 


The  coaxial  injector  was  placed  in  an  inner  chamber  within  the  main  test  chamber.  This  unit  was 
protruding  into  the  chamber  by  5  mm  through  a  17  mm  in  diameter  central  orifice  in  the  top  wall.  This 
inner  chamber  was  created  to  augment  the  level  of  acoustic  modulation  acting  on  the  coaxial  flow.  The 
inner  chamber  had  a  height  of  6.6  cm,  a  width  of  7.6  cm  and  a  depth  of  1.3  cm.  The  flow  was  exhausted 
through  another  orifice  at  the  center  of  the  bottom  wall. 

Injection  temperatures  were  measured  with  unshielded  traveling  type  E  thermocouples  with  a  bead 
diameter  of  0.1  mm.  The  accuracy  of  the  thermocouples  used  in  the  study  was  checked  with  an  RTD 
and  found  to  be  within  1  K.  A  miniature  pressure  transducer  (Kulite  Semiconductor  Products,  models 
CCQ-062-1000A  and  CCQ-093-750A)  was  placed  next  to  the  thermocouple  using  a  small  metallic  post 
for  support,  and  recorded  pressure  fluctuations  at  a  sampling  frequency  of  20  kHz.  These  pressure 
transducers  having  an  absolute  pressure  range  of  either  6.9  MPa  or  5.2  MPa,  respectively  were  used  to 
used  to  scan  the  wavefield  induced  by  the  acoustic  driver  units  prior  to  flow  experiments. 
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The  coaxial  jet  structure  was  observed  with  backlighting  (shadowgraph)  imaging.  The  flow  was  illu¬ 
minated  with  a  light  source  (Newport  variable  power  arc  lamp  set  at  160W  or  300W)  and  the  transmitted 
light  was  recorded  with  a  Phantom  v7.1  high-speed  camera  capable  of  recording  up  to  160,000  frames  per 
second  at  a  resolution  of  32  pixels  by  32  pixels.  The  camera  equipped  with  AF  Nikkor  35-105  mm  lens 
and  a  Nikon  No.  1  Close-Up  lenscan  be  seen  at  the  bottom  center  in  Fig.  3.1(a).  The  image  resolution 
in  the  present  experiments  varied  from  128  pixels  x  224  pixels  to  196  pixels  x  400  pixels,  depending 
on  chamber  pressure  and  outer  to  inner  momentum  flux  ratio.  The  images  were  in  a  grayscale  scale 
with  each  pixel  taking  a  value  between  0  and  255  for  different  intensities  of  gray  with  as  extreme  values 
white  and  black.  Each  pixel  represents  an  approximate  area  of  0.08  mm  x  0.08  mm.  Depending  on  the 
resolution  chosen,  the  framing  rate  was  20,  25  or  41  kHz.  The  exposure  time  was  generally  1-2  fi s  and 
the  number  of  images  saved  per  run  was  1000. 

The  images  recorded  were  then  processed  using  a  Matlab  subroutine  based  on  the  Otsu  technique  [44] . 
The  subroutine  takes  a  grayscale  image  and  finds  a  pixel  threshold  value.  All  the  pixel  values  above  this 
threshold  are  converted  to  black  and  the  remaining  pixels  to  white.  Then  the  subroutine  measures  the 
length  of  the  black  region  extending  from  the  exit  of  the  coaxial  jet.  This  is  referred  to  as  the  “dark  core 
length” .  A  “mean  dark  core  length”  l^c  is  then  obtained  by  averaging  over  the  set  of  images  recorded  by 
the  camera. 

Acoustic  modulation  of  this  setup  was  obtained  with  two  piezo-sirens  placed  at  each  end  of  the 
chamber  at  a  distance  of  34  cm  from  the  jet  axis,  custom-designed  for  the  Air  Force  Research  Laboratory 
by  Hersh  Acoustical  Engineering,  Inc  (the  piezo-sirens  can  be  seen  at  both  sides  of  the  chamber  in  Fig. 
3.1(a)).  A  Fluke  100  MS/s  arbitrary  waveform  generator  (model  292)  was  used  to  produce  two  sinusoidal 
waves  with  the  same  frequency  but  with  a  prescribed  phase  between  them.  The  signals  were  then  sent 
to  two  amplifiers  (Krohn-Hite  model  7500  and  a  Trek  model  PZD2000A),  one  for  each  piezo-siren.  The 
amplified  signals  were  in  the  200  to  540  V  range.  To  accommodate  for  the  rectangular  chamber,  a 
waveguide  with  a  catenary  contour,  see  upper  left  corner  of  Fig.  3.1(a),  was  used  to  transmit  the  waves 
from  the  circular  cross-section  of  each  siren  to  a  rectangular  cross-section.  The  maximum  root-mean- 
square  acoustic  pressure  fluctuations  generated  by  the  two  piezo-ceramic  acoustic  sources  in  the  inner 
chamber,  measured  using  the  three  differential  Kulite  pressure  transducers  (model  XCQ-093-25D)  flush 
mounted  on  a  lateral  chamber  wall,  varied  from  8.0  to  22  kPa  and  were  produced  in  the  2.9  to  3.1 
kHz  range.  These  operating  conditions  were  obtained  by  manually  varying  the  frequency  on  the  signal 
generator  and  finding  the  highest  possible  pressure  amplitude. 

The  phase  difference  between  the  signals  sent  to  the  piezo-siren  elements  was  varied  to  expose  the 
coaxial  jet  to  different  effective  positions  relative  to  the  pressure  node  or  antinode  of  the  acoustic  field. 
When  the  two  piezo-siren  elements  produce  waves  with  a  zero  degree  phase  angle  between  them,  the 
motion  of  the  piezo-siren  radiating  membranes  is  synchronized  and  in  opposite  directions.  This  produces 
conditions  of  high  pressure  perturbations  and  low  velocity  fluctuations  at  the  center  of  the  chamber, 
where  the  coaxial  jet  is  located  which  correspond  to  a  pressure  antinode  or  velocity  node.  In  contrast, 
when  the  two  drivers  present  a  180-degree  phase  difference,  then  the  radiating  membranes  move  in  the 
same  direction,  generating  high  velocity  fluctuations  with  very  small  variations  in  pressure  at  the  jet 
location  producing  a  pressure  node  or  velocity  antinode. 
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Chapter  4 

Balance  equations  and  models 


At  this  stage  it  is  worth  briefly  reviewing  the  spatially  filtered  balance  equations  and  introduce  the  real 
gas  models  and  subgrid  scale  closures  used  in  the  large  eddy  simulations. 


4.1  Governing  equations 


A  single  species  flow  configuration  is  investigated  in  the  present  study.  In  this  situation,  the  vector  of 
conservative  variables  for  a  compressible  flow  is  w  =  (p,  pit,  pE)T ,  where  p  is  the  density,  u  the  velocity 
vector  and  E  the  total  energy  ( E  is  the  sum  of  internal  energy  es  and  kinetic  energy  e&  =  1/2  J2u?). 
The  present  large-eddy  simulations  are  carried  out  by  integrating  the  mass  weighted  spatially  filtered 
Navier-Stokes  equations.  In  the  following  balance  equations  (0  designates  a  spatially  filtered  variable 
while  0  is  a  mass-weighted  (Favre)  spatially  filtered  variable  : 
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In  Eqs.  4. 1-4.3,  p  is  the  pressure,  x  the  spatial  coordinates  vector  and  t  is  the  time.  The  laminar  viscous 
stress  tensor  r  and  heat  flux  vector  q  are  expressed  as  linear  functions  of  the  strain  rates  and  temperature 
gradient  respectively.  Laminar  viscosity  and  heat  conductivity  coefficients  are  determined  with  the  Chung 
et  al.  method  [12].  For  the  single  species  case  investigated  in  the  present  article,  no  additional  diffusion 
terms  are  needed  since  the  Soret  and  Dufour  effects  naturally  vanish.  The  modeling  of  the  subgrid-scale 
stress  tensor  r\-  and  heat  flux  vector  ql-  requires  some  discussion.  It  is  known  from  recent  studies  that 
subgrid  closure  rules  for  real  gases  differ  from  those  used  for  perfect  gas  flows.  In  principle  one  should 
take  into  account  additional  terms  arising  from  the  nonlinearity  of  the  state  equation,  and  from  the  large 
variations  in  transport  properties  as  explained  in  some  recent  studies  by  Bellan  et  al.  [2] ,  Selle  et  al.  [63] 
and  Taskinoglu  and  Bellan  [68].  However  the  modeling  issues  are  not  completely  settled.  It  is  also  known 
from  other  calculations  that  the  low  pressure  subgrid  scale  models  yield  suitable  results  for  transcritical 
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and  supercritical  round  jets  [59,  62].  It  was  decided  for  these  reasons  to  keep  these  low  pressure  closures 
for  the  subgrid  scale  models  as  described  in  more  detail  in  Sec.  4.3. 


4.2  Real  gas  thermodynamics  and  equation  of  state 


Real  gas  thermodynamics  for  high  pressure  flow  calculations  is  envisaged  in  many  recent  studies  where 
accurate  descriptions  of  the  fluid  state  are  derived  from  various  types  of  equations  [1,  71].  The  present 
simulations  rely  on  the  Peng-Robinson  equation  of  state  [45]  (Eq.  4.4)  which  is  used  in  the  analysis  to 
model  the  departure  from  an  ideal-gas  behavior.  This  EOS  is  less  accurate  than  the  modified  version  but 
it  is  also  simpler  and  requires  a  smaller  amount  of  calculations.  The  fluid  state  is  expressed  in  the  form  : 


P  = 


prT 

1  —  pb  1  +  2  pb  —  p2b2 


P2a(T) 


(4.4) 


where  T  is  the  temperature,  r  =  R/W  where  R  is  the  perfect-gas  constant  and  W  the  molar  mass.  The 
Peng-Robinson  coefficients  a(T )  and  b  for  a  single-species  fluid  are  [47]: 
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with 


c  =  0.37464  +  1.54226a;  —  0.26992a; 


(4.7) 


where  Tc  =  126.2  K  and  pc  =  3.396  MPa  are  the  critical  temperature  and  pressure  of  nitrogen,  u)  = 
0.0372  represents  its  acentric  factor  and  TR  =  T /Tc  is  the  reduced  temperature.  This  equation  offers  a 
good  trade-off  between  computational  cost  and  precision.  The  pressure  dependence  of  thermodynamic 
functions  and  coefficients  are  subsequently  derived  from  the  EOS  [33,  38,  71].  Departure  energy  Aes  and 
specific  heat  A cv  at  constant  volume  are  given  by  [20,  47]: 


where  Aes  ■=  es(p,  T)  —  e^(T)  and  A cv  =  cv(p,T)  —  c®(T).  The  enthalpy  and  specific  heat  at  constant 
pressure  are  then  deduced  from  : 
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The  thermal  expansion  coefficient  a  and  the  isothermal  compressibility  coefficient  /3  appearing  in  the 
previous  expressions  are  given  by: 
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The  low-pressure  references  {i.e.  and  c£)  are  obtained  from  the  JANAF  Thermodynamical  Tables  [7], 
and  only  depend  on  temperature. 
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4.3  Subgrid-scale  models  and  assumptions 


As  indicated  previously,  the  present  simulations  rely  on  a  standard  (low  pressure)  subgrid  scale  closure, 
the  WALE  model,  in  which  the  stress  tensor  rj-  is  conveniently  expressed  in  terms  of  the  spatially  resolved 

strain  rates  Sij  [35]  : 
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In  the  previous  expression  the  WALE  model  constant  Cw  =  0.4929  and  gij  is  the  irrotational  part  of  the 
stress  tensor  : 


—  UUji  .  ~  ~  /  a  tr\ 

9ij  =  7^77  and  g{j  =  gikgkj  (4.15) 

The  WALE  model  distinguishes  velocity  gradients  associated  with  the  rotational  and  pure  shear  compo¬ 
nents  of  the  velocity  gradients  in  order  to  distinguish  regions  where  turbulence  is  transitional  and  regions 
where  it  is  is  fully  developed.  It  is  well  suited  to  the  treatment  of  shear  flows  of  the  type  considered  in 
this  article. 

The  subgrid  scale  heat  flux  is  modeled  using  a  standard  gradient  transport  assumption  : 

=  -A«g  (4.16, 

The  subgrid  scale  conductivity  is  cast  in  the  form  Xt  =  cpgt/ Prt  where  fit  =  pvt  and  Prt  respectively 
designate  the  subgrid  dynamic  viscosity  and  turbulent  Prandtl  number.  As  in  low  pressure  simulations 
a  constant  value  Pr*  =  0.7  is  used  in  the  present  calculations. 
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Chapter  5 

Numerical  aspects 


The  mass  weighted  spatially  filtered  compressible  Navier- Stokes  equations  in  combination  with  the  subgrid 
scale  models  described  previously  are  integrated  in  the  AVBP  flow  solver  [60,  34].  This  code  has  already 
been  used  to  calculate  low  and  moderate  pressure  reactive  flows  of  multicomponent  mixtures  of  perfect 
gases.  The  spatial  discretization  on  unstructured  or  hybrid  meshes  facilitates  applications  in  which  the 
geometry  is  complex,  a  feature  which  is  quite  useful  if  one  wishes  to  analyze  practical  configurations. 
The  LES  of  turbulent  flows  requires  a  low-dissipation  algorithm  [21,  18]  which  is  often  obtained  with 
centered  discretization  schemes.  The  present  integration  method  relies  on  a  Taylor- Galer kin  weighted 
residual  central  distribution  scheme,  called  TTG4A.  This  scheme  is  third-order  in  space  and  fourth-order 
in  time  [50]. 

The  real  gas  equation  of  state,  thermodynamics  and  transport  coefficients  have  been  implemented  in 
AVBP.  The  convective  fluxes  Jacobian  matrices,  used  by  the  scheme,  are  expressed  in  terms  of  real  gas 
thermodynamics  to  preserve  the  overall  consistency  of  the  code.  A  fully  consistent  treatment  of  associ¬ 
ated  boundary  conditions  is  based  on  the  characteristic  wave  decomposition  method  NSCBC  [46,  34]  and 
includes  specific  expressions  provided  for  real  gas  thermodynamics  in  [39].  The  highly  nonlinear  thermo¬ 
dynamics  of  the  transcritical  fluid  stream  induces  large  density  gradients  between  the  dense  transcritical 
fluid  and  the  surrounding  gaseous  stream  which  require  a  specific  stabilization  procedure.  This  relies  on 
artificial  viscosity  and  is  used  when  unresolved  gradients  are  detected  following  a  method  described  in 
[59].  The  real  gas  AVBP  flow  solver  has  already  been  used  to  simulate  transcritical  cryogenic  flames  of 
various  types  [48,  58,  57]. 


5.1  Computational  domain  and  mesh 


The  geometry  of  the  computational  domain  essentially  matches  the  experimental  configuration.  A  three 
dimensional  view  of  the  external  boundaries  is  given  in  Fig.  5.1(a).  The  longitudinal  section  of  the 
chamber  is  rectangular  (59.4x76.0  mm)  and  its  width  is  12.7  mm.  The  fluid  streams  are  injected  by  a 
coaxial  element  which  protrudes  by  5  mm  into  the  chamber  (Fig.  5.1b).  This  injector  is  located  in  the 
center  of  a  17  mm-diameter  hole  in  the  inlet  wall.  No  co-flow  is  directly  injected  through  this  orifice 
during  the  experiment  but  nitrogen  gas  may  traverse  this  aperture  as  a  result  of  entrainment  by  the 
coaxial  streams  giving  rise  to  a  finite  flow  velocity  in  this  region.  The  entrained  fluid  originates  from 
a  large  reservoir  defined  on  the  upstream  side  of  the  central  hole.  Details  about  the  injector  geometry 
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are  given  in  Fig.  5.2(a).  The  inner  injector  diameter,  inner  duct  thickness,  annular  duct  thickness  and 
outer  duct  diameter  are  respectively  di  =  0.51  mm,  li  =  0.54  mm,  he  =  0.415  mm  and  dQ  =  2.42  mm. 
The  mesh  comprises  2  100  000  nodes  corresponding  to  10  000  000  tetrahedra.  It  is  highly  refined  near 


59.4  mm 


Figure  5.1:  (a)  Three-dimensional  visualization  of  the  reservoir,  (b)  Longitudinal  cut  of  the  mesh. 


the  injector  (Fig.  5.2(b)),  with  a  constant  characteristic  size  of  0.032  mm  on  a  distance  of  10  inner  jet 
diameters.  The  mesh  is  then  slowly  coarsened  towards  the  domain  outlet  (Fig.  5.1(b)). 


Figure  5.2:  (a)  Closer  view  of  the  injector,  (b)  Longitudinal  cut  of  the  mesh  (35  inner  injector  diameters). 


5.2  Injection  characteristics  and  boundary  conditions 


Three  test  conditions  are  considered  in  this  analysis  designated  as  “N2”,“N6”  and  “N8”.  Case  “N6” 
serves  as  a  reference  simulation.  The  two  other  cases,  designated  as  “N2”  and  “N8”,  are  also  calculated 
in  order  to  analyze  effects  of  the  momentum  flux  ratio  on  the  flow  behavior  (Ch.  6).  The  three  cases 
differ  by  their  injection  temperatures  and  velocities,  yielding  momentum  flux  ratios  J  of  1,  2.6  and  9.1 
for  cases  “N2”,  “N6”  and  “N8”  respectively.  The  corresponding  density  ratios  S  are  5.84,  5.74  and  3.22 
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respectively.  Injection  conditions  are  gathered  in  Tab.  3.1  and  Tab.  3.2  which  respectively  correspond  to 
the  outer  and  inner  streams.  Injection  conditions  are  placed  on  thermodynamic  plots  on  Fig.  5.3.  The 
density  and  volumetric  internal  energy  notably  depart  from  a  perfect  gas  behavior  when  temperature  is 
close  and  below  the  critical  value.  The  influence  of  acoustic  perturbations  will  only  be  investigated  in  a 
single  case  corresponding  to“N6”  (Ch.  7). 

In  all  cases,  pressure  is  maintained  at  the  outlet  boundary  using  non-reflecting  characteristic  conditions 
at  an  adapted  level  of  3.56  MPa  and  the  chamber  and  reservoir  are  initially  filled  with  nitrogen  at  the 
same  pressure  and  a  temperature  213  K  (the  corresponding  density  is  59.6  kg  m~ 3).  All  solid  boundaries 
are  treated  as  adiabatic  slip  walls.  Velocity  perturbations  (with  an  amplitude  of  3  %  of  the  mean  flow) 
are  added  at  injection  in  order  to  represent  effects  of  turbulent  fluctuations  using  a  procedure  described 
in  [6,  66,  17]. 


Figure  5.3:  Thermodynamic  state  calculated  with  the  PR  EOS  ( — )  and  the  perfect  gas  EOS  (•  •  •).  (a) 
Density  in  terms  of  temperature,  (b)  Volumetric  internal  energy  (i.e.  pes)  in  terms  of  temperature.  The 
pressure  is  constant  and  equal  to  3.56  MPa.  Symbols:  o  inner  jet  conditions;  x  outer  stream  conditions; 
•  chamber  conditions. 
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N2 

1.05 

0.42 

5.86 

2.48 

N6 

3.05 

0.72 

5.65 

4.15 

N8 

9.3 

1.70 

3.22 

5.48 

Table  5.1:  Injection  characteristics  of  the  simulated  cases.  J  is  the  momentum  flux  ratio  between  the 
outer  and  the  inner  jets,  M  the  mass  flux  ratio,  S  the  density  ratio  and  V  is  the  velocity  ratio,  pi,  pe,  iq 
and  ue  are  the  inner  and  outer  jet  densities  and  velocities. 


5.3  Simulation  procedure 


Simulations  are  carried  out  on  a  coarse  mesh,  with  an  initial  solution  corresponding  to  a  uniform  fluid 
at  rest.  The  flow  is  then  qualitatively  established  after  three  convective  periods.  Results  are  then 
interpolated  on  the  final  mesh.  Data  averaging  can  be  started  after  one  transient  period.  The  simulated 
physical  time  is  given  in  Tab.  5.2  together  with  the  number  of  simulated  convective  periods.  Finally,  the 
CPU  time  (equal  to  the  product  ncore  x  rihours )  needed  to  obtain  the  average  solution  (which  does  not 
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include  the  CPU  time  needed  for  the  transient  flow  and  the  simulation  on  the  coarse  mesh)  is  also  given 
in  this  table. 


Case 

A ta  [ms] 

A ta/ru t 

A  ta/rUe 

CPU  time  [h] 

N2 

14.3 

1.9 

4.8 

35  000 

N6 

22.7 

4.4 

18.2 

45  000 

N8 

10.4 

4.0 

21.2 

20  000 

Table  5.2:  Averaging  time  for  cases  “N2”,  “N6”  and  “N8”.  A ta  is  the  physical  averaging  time.  ru. 
represents  one  convective  time  over  35  di  for  the  inner  jet  and  rUe  represents  one  convective  time  over  35 
di  for  the  outer  jet.  rUi  =  (35 di)/ui  ;  rUe  =  (35 di)/ue. 
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Chapter  6 

Effect  of  momentum  flux  ratio 


It  is  first  natural  to  consider  the  influence  of  the  momentum  flux  ratio  on  the  system  behavior.  Flow 
patterns  calculated  numerically  and  obtained  experimentally  are  first  compared  in  Sec.  6.1.  This  is 
followed  by  a  description  of  the  mean  flow  in  Sec.  6.2  and  by  an  analysis  of  flow  dynamics  in  Sec.  6.3. 


6.1  Instantaneous  flow  patterns 


Instantaneous  backlighting  images  obtained  experimentally  are  shown  in  Fig.  6.1(a),(c),(e).  In  the  dark 
core  light  is  deviated  by  changes  in  the  refraction  index.  This  in  turn  delineates  regions  where  the  density 
takes  large  values  (the  deviation  can  be  roughly  linked  to  the  second  derivative  of  the  density  integrated 
on  the  line  of  sight).  It  is  generally  admitted  that  the  dark  region  represents  the  position  of  the  dense 
core,  and  qualitatively  highlights  its  external  boundaries  giving  an  insight  on  the  main  structures  present 
in  the  flow  [43].  As  the  momentum  flux  ratio  is  increased,  the  inner  jet  axial  length  is  reduced.  Case 
“N2”  exhibits  large  scale  structures  over  about  10  inner  injector  diameters  downstream  of  the  injector. 
These  structures  are  not  identifiable  in  cases  “N6”  and  “N8”.  The  latter  is  characterized  by  a  sudden 
termination  of  the  inner  jet  at  a  few  inner  jet  diameters  from  the  exit  plane.  In  contrast  with  cases  “N2” 
and  “N6”,  the  inner  stream  appears  to  be  quickly  mixed  with  the  surrounding  fluid  and  the  central  jet 
can  no  longer  be  clearly  identified  downstream  of  this  abrupt  termination. 

Longitudinal  slices  of  instantaneous  density  distributions  extracted  from  simulations  are  shown  in 
Fig.  6.1(b),(d),(f)  which  respectively  correspond  to  “N2”,  “N6”  and  “N8”  cases.  Comparisons  between 
these  figures  and  the  experimental  backlighting  images  are  admittedly  qualitative.  Nevertheless,  general 
experimental  features  corresponding  to  these  different  cases  are  replicated.  One  finds  that  case  “N2” 
features  large  scale  structures  which  are  not  present  in  cases  “N6”  and  “N8”.  The  abrupt  termination 
noticed  in  case  “N8”  is  retrieved  numerically.  This  is  also  made  more  apparent  in  the  density  iso¬ 
surfaces  corresponding  to  p  —  0.1  {pi  —  pe)  +  pe)  plotted  in  Fig.  6.2.  The  large  scale  structures  observed 
experimentally  characterize  case  “N2”,  and  disappear  as  the  momentum  flux  ratio  is  increased  (cases 
“N6”  and  “N8”). 
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□ 

(f) 


Figure  6.1:  Comparison  between  experimental  visualizations  (backlighting)  and  typical  instantaneous 
density  distributions  obtained  from  simulations  (white:  60  kg  m-3  ;  black:  inner  jet  injection  density  (Tab. 
3.2)  ;  logarithmic  scale).  Case  “N2”  (J=1.05)  :  (a)  Experimental  backlighting  image,  (b)  Calculated 
density.  Case  “N6”  (J=3.05)  :  (c)  Experimental  backlighting  image,  (d)  Calculated  density.  Case  “ N8 ” 
(J=9.3)  :  (e)  Experimental  backlighting  image,  (f)^galculated  density. 


□ 


J 


m 


□ 

(b) 

□ 


□ 

(c) 

Figure  6.2:  LES  results.  Density  iso-surface  corresponding  to  p  =  0.1  {pi  —  pe)  +  pe  colored  by  the  axial 
velocity,  (a)  case  “N2”  (J=1.05)  umin=- 2  m  s-2,  ^ma^=6.5  m  s-2,  (b)  case  “N6”  (J=3.05)  umin=-4:  m 
s-2,  urnax= 16  m  s-2,  (c)  case  “N8”  (J=9.3)  umin=-5  m  s-2,  umax=40  m  s-2.  The  plot  corresponds  to 
20  inner  injector  diameters. 
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6.2  Mean  flow  description 


It  is  first  interesting  to  examine  effects  of  the  momentum  flux  ratio  on  the  mean  “dark  core”  length.  Since 
the  density  field  cannot  be  directly  compared  to  the  experimental  light  distribution  images,  it  is  natural 
to  define  a  set  of  density  length  scales  and  compare  them  with  the  dark  core  length  deduced  from  the 
data.  The  jet  core  characteristic  lengths  are  defined  by  considering  the  axial  density  distribution  : 


;0.99 

V 


-  -0.99 

-  xp 


- 


=  x"*-x  o 


iy=x^-x  0 


(6.1) 


where  xo  is  the  axial  position  of  the  inner  injector  exit,  x°p 99 ,  x°p 5  and  x°pA  are  positions  on  the  jet  axis 
where  density  is  respectively  equal  to  99  %  of  its  initial  value,  to  0.5 (pi  —  pe)  +  pe  or  to  0.1  {pi  —  pe)  +  pe. 
These  three  lengths  are  shown  in  Fig.  6.3(a)  together  with  mean  density  profiles  on  the  centerline. 


It  is  also  useful  to  define  at  this  point  another  length  scale  based  on  the  axial  mean  velocity  profile. 
The  central  jets  all  feature  an  initial  decrease  in  the  axial  velocity  on  the  centerline  as  can  be  seen  in  Fig. 
6.3(b).  The  minimum  value  is  reached  at  xu  which  is  used  to  define  a  new  length  scale  lu  =  xu  —  xq. 


Characteristic  lengths  are  gathered  in  Tab.  6.1  for  comparison  with  experimental  measurements  of 
the  “dark  core”  length.  The  latter  values  are  deduced  from  backlighting  images  by  measuring  the  “dark 
region”  size  from  average  images.  Again,  comparison  between  experiment  and  simulation  is  limited.  As 
already  noted  in  [43],  lengths  deduced  from  experimental  backlighting  data  are  generally  larger  than 
those  deduced  from  density  measurements.  It  is  however  clear  that  there  are  similarities  between  trends 
in  simulations  and  experiments.  In  agreement  with  experiments,  numerical  results  indicate  a  reduction 
in  the  density  length  scales  as  the  momentum  flux  ratio  is  increased.  The  velocity  length  scale  is  reduced 
and  the  amplitude  of  the  velocity  defect  is  enhanced  as  the  momentum  flux  ratio  is  increased.  The  flow 


Figure  6.3:  Centerline  profiles  deduced  by  averaging  LES  results,  (a)  Density,  (b)  Axial  velocity  normal¬ 
ized  by  the  outer  injection  velocity  u *  =  u/ue.  -  :  Case  “N2”,  —  •  —  :  Case  “N6”  and - :  Case 

“N8”. 


may  also  be  characterized  by  the  distributions  of  mean  axial  velocity  plotted  in  Fig.  6.4.  The  outer  flow 
contracts  towards  the  central  stream  and  pinches  the  central  jet  to  a  point  where  this  stream  disappears. 
In  the  three  cases,  a  small  back- flow  is  generated  behind  the  injector  lips.  As  the  momentum  flux  ratio 
is  increased  beyond  a  certain  critical  value  (case  “N8”),  the  central  jet  is  surrounded  by  a  negative  axial 
velocity  region  (Fig.  6.4c). 

Case  “N8”  exhibits  a  large  negative  value,  as  observed  previously  in  Fig.  6.4(c).  Further  downstream, 
the  flows  behave  in  a  self-similar  fashion  in  terms  of  axial  velocity  normalized  by  the  outer  jet  injection 
velocity.  Beyond  this  section,  the  impact  of  the  inner  jet  on  the  flow  disappears  as  this  jet  is  mixed  with 
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Case 

Ca9M 

C6M 

CpA/di 

IJ  di 

ldc/di 

N2 

2.57 

4.40 

8.66 

3.8 

15.5 

N6 

1.81 

3.37 

5.88 

2.95 

6.0 

N8 

0.85 

1.17 

1.67 

1.96 

2.5 

Table  6.1:  Characteristic  lengths  deduced  from  simulations  based  on  axial  density  and  velocity  profiles 
and  from  experimental  measurement  of  the  dark  core  size. 


(a) 


(b) 


(c) 

Figure  6.4:  Distribution  of  axial  velocity  (white:  minimum  ;  black:  maximum).  This  map  corresponds  to 
35  inner  injector  diameters.  —  iso-contour  of  zero  axial  velocity,  (a)  case  “ N2 ”  ( J=1.05 ),  (b)  case  “N6” 
(J =3. 05),  (c)  case  “N8 ”  (J=9.3). 
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the  outer  stream.  Since  the  mass  flow  rate  of  the  outer  jet  is  much  higher  than  that  of  the  inner  jet,  the 
outer  jet  dominates  the  flow.  Interestingly,  the  three  profiles  merge  at  x  >  15 indicating  a  self-similar 
behavior  when  profiles  are  plotted  in  terms  of  normalized  variables. 

Radial  profiles  of  normalized  density  and  axial  velocity  are  plotted  in  Fig.  6.5  for  case  “N6”  (case 
“N2”  is  not  shown  here  because  its  behavior  is  quite  similar).  The  first  profiles  lie  in  the  potential  core 
region  (x  <  x°p "),  and  feature  constant  levels  of  density  and  velocity  in  that  region.  Negative  values 
of  axial  velocity  indicate  that  back- flow  is  generated  behind  the  lips.  The  density  profile  behind  the 
inner  lip  shows  a  rapid  variation  towards  the  outer  stream.  This  is  attributed  to  an  efficient  mixing  of 
inner  and  outer  streams  associated  with  the  formation  of  vortical  structure  behind  the  inner  lip  and  the 
generation  of  a  back  flow.  In  this  region,  mass  is  extracted  from  the  central  jet  by  vortical  structures  and 
recirculates  upstream.  Dynamics  in  this  region  is  discussed  in  Sec.  6.3.  The  profile  at  3 di  corresponds 
to  x  =  xu.  The  density  stratification  is  still  apparent.  At  x  =  15 di,  velocity  and  density  are  maximum 
on  the  centerline,  and  a  jet-like  profile  is  retrieved.  Further  downstream,  the  flow  behaves  like  a  single 
round  jet  and  develops  in  a  self-similar  fashion.  The  half  width  at  half  maximum  (HWHM)  of  density 
Lp  is  plotted  in  Fig.  6.6(a).  For  x  >  20 di,  this  quantity  evolves  quasi- linearly  Lp  =  Ejf  (x  —  xv),  where 
E ^  is  the  spreading  rate  of  the  jet  and  xv  is  a  virtual  origin.  Linear  regressions  for  20 di  <x<  35 di  give 
spreading  rates  of  0.12,  0.11  and  0.11  for  cases  “N2”,  “N6”  and  “N8”,  respectively.  These  values  are  in 
agreement  with  those  given  for  example  by  Chen  and  Rodi  [10]  for  a  single  jet  and  with  results  obtained 
from  previous  high  pressure  simulations  of  single  species  round  jets  [59,  76].  In  this  region  of  the  flow 
the  density  profile  takes  a  gaussian  shape  (Fig.  6.6b).  A  linear  evolution  of  the  density  decay  is  also 
observed,  but  not  shown  here,  confirming  that  the  flow  is  self-similar  in  the  downstream  region.  Since 


Figure  6.5:  LES  results ,  Case  “N6”.  (a)  Radial  profiles  of  normalized  density  p*  =  (p  —  poo) /{pc  ~  Poo) 
(b)  Radial  profiles  of  axial  velocity. 


the  flow  develops  in  a  confined  domain  one  may  wonder  if  it  is  influenced  by  recirculation  in  the  chamber. 
Examination  of  velocity  distributions  in  the  domain  indicates  that  the  outer  recirculation  velocities  are 
low  and  that  the  reverse  flow  in  the  outer  region  has  a  limited  influence  on  the  coaxial  jet  structure. 


6.3  Near  injector  flow  dynamics 


Large  eddy  simulations  allow  detailed  analysis  of  the  flow  dynamics.  This  can  be  accomplished  by 
estimating  the  spectral  content  of  fluctuations  in  various  regions  of  the  flow  to  identify  the  preferred 
modes  of  the  system.  To  this  end,  it  is  convenient  to  place  a  set  of  sensors  in  the  flow  and  examine  the 
corresponding  velocity  and  density  time  series.  This  is  illustrated  in  Fig.  6.7.  Since  the  general  flow 
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Figure  6.6:  LES  results  (a)  Half  Width  Half  Max  of  density - Case  “N2”,  —  •  —  Case  “N6”  and - 

Case  “N8”.  (b)  Radial  profiles  of  normalized  density  p'  =  {p  —  poo) /{pc  ~  Poo)- 


structures  of  cases  “N2”  and  “N6”  differ  from  that  corresponding  to  “N8”,  it  is  natural  to  study  these 
two  cases  separately  (Sec.  6.3.1).  Case  “N8”  is  then  considered  in  Sec.  6.3.2. 


Figure  6.7:  The  three  sensor  locations.  The  numerical  probes  are  used  to  measure  characteristic  frequen¬ 
cies  in  the  coaxial  flow,  (a)  Case  “N6”  (J=3.05).  Density  distribution  (white:  60  kg  m-3  ;  black:  410  kg 
m-3  ;  logarithmic  scale).  Sensor  B  is  placed  at  the  end  of  the  inner  jet  potential  core  which  corresponds 
to  different  location  in  cases  “N2”  and  “N6”.  (b)  Case  “N8”  (J  =  9.3).  Density  distribution  (white:  60 
kg  m-3  ;  black:  maximum  ;  logarithmic  scale).  Sensor  C  is  located  in  the  recirculation  region. 


6.3.1  Cases  “N2”  and  “N6” 

To  examine  the  spectral  information  it  is  natural  to  define  a  set  of  Strouhal  numbers  associated  with 
characteristic  scales  and  velocities  in  the  region  of  interest.  The  meaningful  length  scales  in  this  study  are 
the  inner  injector  diameter  di,  the  inner  lip  thickness  li  and  the  annular  jet  thickness  he.  Two  velocities 
characterize  the  system  :  the  inner  and  outer  jet  injection  velocities  (iq  and  ue ,  respectively).  The 
resulting  Strouhal  numbers  and  the  probes  used  to  measure  frequencies  are  shown  in  Tab.  6.2.  Probes 
A  and  B  are  used  to  evaluate  the  fundamental  frequencies  of  the  outer  and  inner  streams.  Probe  C 
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provides  the  dominant  frequency  associated  with  eddy  structures  released  behind  the  injector  lip  and 
associated  with  the  outer  jet.  From  instantaneous  snapshots  (see  Fig.  6.7  for  example)  and  animations, 
these  structures  are  found  to  originate  from  the  lip  corner  on  the  outer-jet  side.  Positions  of  the  different 
probes  are  shown  in  Fig.  6.7  (a).  Note  that  the  axial  location  of  probe  B  is  not  the  same  in  cases  “N2” 
and  “N6”,  because  this  probe  is  placed  at  the  end  of  the  potential  core.  Power  spectral  densities  are 
estimated  using  Welch’s  method  of  periodograms.  Statistical  averaging  is  obtained  by  segmenting  the 
finite  data  set  in  M  overlapping  blocks  (a  50  %  overlap  is  used  to  augment  the  number  of  periodograms 
used  in  the  averaging  procedure).  Because  the  number  of  samples  is  limited  the  averaging  is  carried 
out  over  M  =  8  blocks  in  the  outer  flow  and  on  M  =  4  blocks  in  the  inner  jet  (power  is  concentrated 
in  the  low  frequency  range  and  it  is  important  to  augment  the  resolution  to  distinguish  the  maximum 
frequency).  The  frequency  resolution  corresponding  to  the  different  sensors  is  given  in  Tab.  6.2. 


St  =  ( fl)/u 

l 

u 

Probe  (Fig.  6.7(a)) 

A  /  [Hz] 

Ste 

he  (0.41  mm) 

ue 

A 

N2:  375  ;  N6:  250 

St 1 

di  (0.51  mm) 

B 

N2:  187  ;  N6:  125 

St1 

li  (0.54  mm) 

Ue 

C 

N2:  375  ;  N6:  250 

Table  6.2:  Definition  of  characteristic  Strouhal  numbers.  The  characteristic  length  and  velocity  l  and  u 
are  given  in  the  second  and  third  columns.  The  frequency  /  corresponds  to  the  power  spectral  density 
maximum  estimated  from  the  simulations.  A /  is  the  frequency  resolution. 


The  outer  and  central  streams  are  destabilized  at  a  few  diameters  di  from  the  injector  exit.  The 
measured  frequency  for  the  outer  jet  (probe  A)  is  3  600  Hz  for  case  “N2”  and  8  550  Hz  for  case  “N6” 
(Fig.  6.8),  which  correspond  to  St^2=0.25  and  St^6=0.25.  The  central  high-density  jet  presents  a 
fundamental  instability  at  1  620  Hz  and  1  740  Hz  for  cases  “N2”  and  “N6”  (Fig.  6.9)  respectively 
(probe  B).  The  associated  Strouhal  numbers  are  then  St]V2=0.34  and  St]V6=0.26.  For  both  cases,  the 
fundamental  mode  of  the  inner  and  outer  streams  fall  in  the  typical  range  for  round  jets  (0.2  -  0.4). 

Outer  jet 


Figure  6.8:  Power  spectral  density  of  the  radial  velocity  in  the  outer  jet  (probe  A).  Averaging  is  carried 
out  on  M  =  8  blocks  with  a  50%  overlap,  (a)  Cases  “N2”  ( J=1.05),  (b)  Case  “N6”  (J=3.05). 


The  mixing  layer  between  the  two  jets  features  eddy  structures  generated  downstream  of  the  central 
lip.  Their  dynamics  is  similar  to  that  observed  for  vortices  shed  behind  a  backward- facing  step.  The 
characteristic  frequencies  are  found  to  be  1  620  Hz  and  3600  Hz,  for  cases  “N2”  and  “N6”  (Fig.  6.10), 
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Inner  jet 


Figure  6.9:  Power  spectral  density  of  radial  velocity  in  the  inner  jet  (probe  B).  Averaging  is  carried  out 
on  M  =  4  blocks  with  a  50%  overlap,  (a)  Cases  “N2”  (J=1.05),  (b)  Case  “N6”  (J=3.05). 


respectively.  Such  frequencies  correspond  to  Strouhal  numbers  St^— 0.15  and  St5y6=0.23  for  cases  “N2” 
and  “N6”  which  are  slightly  higher  than  the  commonly  determined  values  for  backward-facing  step 
instabilities  (St=O(0.1)  [70]).  This  could  be  attributed  to  a  reduction  of  the  “effective”  lip  thickness 
as  shear  is  increased.  It  turns  out  that  the  backward-facing  step  vortical  structures  induce  a  back-flow 
behind  the  lip  (Fig.  6.4).  High  density  fluid  from  the  central  jet  is  entrained  and  recirculates  behind  the 
lip.  One  may  then  consider  that  this  results  in  a  reduced  “effective”  lip  thickness  which  intervenes  in  the 
formation  of  the  vortical  structures,  thus  limiting  their  spatial  development  and  increasing  the  instability 
frequency.  All  the  Strouhal  numbers  obtained  are  gathered  in  Tab.  6.3.  The  frequency  obtained  for  case 
“N2”  is  interestingly  close  to  that  associated  with  the  inner  jet  preferred  mode.  This  observation  may 
explain  why  large  scale  structures  are  exclusively  observed  in  this  case  (Fig.  6.1). 


Inner  lip 


Figure  6.10:  Power  spectral  density  of  radial  velocity  in  the  mixing  layer  behind  the  inner  lip  (probe  C). 
Averaging  is  carried  out  over  M  =  8  blocks  with  a  50%  overlap,  (a)  Cases  “N2”  (J=1.05),  (b)  Case  “N6” 
(J=3.05). 
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Case 

Ste  (probe  A) 

St1  (probe  B) 

St1  (probe  C) 

N2 

0.25 

0.34 

0.15 

N6 

0.25 

0.26 

0.23 

Table  6.3:  Characteristic  Strouhal  numbers  obtained  from  the  simulations  for  cases  “N2”  and  “N6”. 

6.3.2  Case  UN8”:  recirculation  regime 

The  mechanism  leading  to  mixing  of  the  high  density  jet  with  the  outer  stream  is  somewhat  different 
for  case  “N8”.  It  is  already  known  that  the  dense  jet  ends  abruptly  and  that  a  recirculation  bubble  is 
established  beyond  that  point.  As  done  previously  for  cases  “N2”  and  “N6”,  characteristic  length  scales 
used  to  define  the  relevant  Strouhal  numbers  are  given  in  Tab.  6.4  and  the  probes  used  to  characterize 
this  case  are  shown  in  Fig.  6.7(b).  Probe  C  is  inside  the  recirculation  region  and  is  located  at  x  =  xu. 
This  case  is  highly  turbulent  and  features  a  broadband  spectrum  so  that  temporal  data  are  arranged 
in  M  =  16  blocks  to  ensure  a  reasonable  statistical  convergence.  A  50  %  overlap  is  used  to  obtain  the 
maximum  possible  averaging  number.  The  resulting  frequency  resolution  is  in  this  case  1  000  Hz.  This 
value  may  seem  quite  high,  but  should  be  compared  to  the  spectral  range  which  extends  to  25  000  Hz. 


St  =  (. fl)/u 

l 

u 

probe  (Fig.  6.7(b)) 

SteN8 

he  (0.41  mm) 

ue 

A 

StlN8 

li  (0.54  mm) 

Ue 

B 

StbN8 

di  +  li  (0.51+0.54  mm) 

Ue 

C 

Table  6.4:  Strouhal  numbers  definitions.  I  and  u  are  the  characteristic  length  and  velocity  used  to  define 
the  Strouhal  numbers.  /  represents  the  frequency  obtained  from  the  simulations  and  A /  is  the  frequency 
resolution. 


The  dominant  frequency  for  the  outer  jet  (probe  A)  is  21  600  Hz  (Fig.  6.11),  which  corresponds  to 
St^g=0.25,  as  also  obtained  for  the  two  other  cases.  The  inner  jet  is  notably  affected  by  the  recirculation 
bubble.  Its  longitudinal  development  is  limited  and  mass  is  ejected  in  the  radial  direction  in  the  form 
of  thin  layers  appearing  in  Fig.  6.7(b).  The  resulting  structure  is  similar  to  a  counterflow.  This  motion 
leads  to  the  formation  of  a  recirculation  ring  of  high  density  fluid  behind  the  inner  lip  with  a  characteristic 
size  of  the  order  of  the  lip  thickness.  Vortex  roll-up  is  observed  behind  the  inner  lip,  in  the  mixing  layer 
established  between  the  inner  and  outer  jets.  These  structures  are  similar  to  those  observed  in  cases 
“N2”  and  “N6”,  but  they  strongly  interact  with  higher  density  layers  emerging  from  the  central  jet  and 
their  spatial  development  is  limited  by  the  recirculation  region.  The  corresponding  Strouhal  number  is 
StlN 8=0.26  (Fig.  6.12(a)).  These  eddy  structures  lead  to  the  formation  of  the  recirculation  region.  The 
measured  frequency  corresponds  to  StbN8= 0.18  (Fig.  6.12(b)).  The  dynamics  of  this  bubble  is  similar  to 
that  observed  behind  a  bluff  body  terminated  by  a  flat  face,  for  which  the  characteristic  Strouhal  number 
is  close  to  0.2.  The  Strouhal  numbers  obtained  in  this  analysis  are  gathered  in  Tab.  6.5. 


St% 8  (probe  A) 

StlN 8  (probe  B) 

StbN 8  (probe  C) 

0.25 

0.26 

0.18 

Table  6.5:  Characteristic  Strouhal  numbers  obtained  from  the  simulations  for  case  “N8”. 
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Outer  jet 


Figure  6.11:  Case  “N8”  (J= 9.3).  Power  spectral  density  of  pressure  on  Probe  A.  Averaging  is  carried 
out  over  M  =  16  blocks  with  50%  overlap. 


Probe  B 


Recirculation  bubble 


Figure  6.12:  Case  “N8”  (J=9.3).  (a)  Power  spectral  density  of  pressure  on  Probe  B,  (b)  Power  spectral 
density  of  pressure  in  the  recirculation  bubble.  Averaging  is  carried  out  over  M  —  16  blocks  with  50% 
overlap. 
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Chapter  7 

Effects  of  transverse  acoustic 
modulation 


Effects  of  an  acoustic  modulation  are  now  investigated  in  case  “N6”.  This  operating  point  is  chosen 
because  the  core  is  sufficient  long  to  allow  an  examination  of  external  modulation  effects.  The  method¬ 
ology  used  to  impose  an  acoustic  excitation  is  described  in  Sec.  7.1.  Results  are  then  compared  with 
experimental  observations  in  Sec.  7.2  and  discussed  in  Sec.  7.3. 


7.1  Modulation  methodology 


An  acoustic  modulation  is  obtained  by  imposing  harmonic  modulations  of  the  normal  velocity  at  the 
outer  boundaries  using  a  method  devised  by  Rey  et  al.  [52]  (Fig.  7.1).  The  modulation  has  the  form  : 

Vi  =  sin(27r/£)  V2  =  sin(27r  ft  +  0)  (7.1) 

where  v\  and  V2  are  the  normal  velocities  at  the  lateral  boundaries,  /  is  the  modulation  frequency  set 
equal  to  3  000  Hz,  0  designates  the  phase  between  the  imposed  modulations.  All  other  boundaries  are 
treated  as  in  the  previous  cases  (in  the  absence  of  modulation).  This  technique  is  first  validated  by 
simulating  an  acoustic  perturbation  in  the  absence  of  the  central  flow.  Longitudinal  distributions  of  the 
rms  pressure  fluctuations  are  shown  in  Fig.  7.2.  When  0  =  0  (Fig.  7.2(a)),  the  normal  velocities  at  the 
boundaries  are  in  phase  and  this  yields  a  maximum  pressure  fluctuation  on  the  injector  axis.  A  pressure 
probe  placed  at  5  mm  from  the  injector  exit  confirms  this  behavior.  As  expected  pressure  fluctuations 
are  observed  with  a  2  %  peak-to-peak  oscillation  with  respect  to  the  mean  chamber  pressure  while  the 
transverse  velocity  remains  unperturbed.  When  the  normal  velocities  at  the  boundaries  are  in  phase 
opposition  (0  =  7 r),  the  pressure  vanishes  on  the  centerline  and  the  transverse  velocity  oscillation  reaches 
its  maximum  (Fig.  7.2(b).  Spectral  analysis  also  indicates  that  the  perturbed  field  is  dominated  by  the 
imposed  modulation  frequency  and  that  the  harmonic  content  is  negligible. 

The  simulated  cases  are  defined  in  Tab.  7.1.  Injection  conditions  correspond  to  case  “N6” .  Modulation 
amplitudes  of  1  and  2%  of  the  mean  chamber  pressure  are  imposed.  These  two  cases  are  studied  for  0  =  0 
and  0  =  7r  respectively  corresponding  to  a  pressure  antinode  and  a  pressure  node  on  the  central  axis. 
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Figure  7.1:  Longitudinal  slice  in  the  computational  domain  showing  the  methodology  retained  to  impose 
acoustic  perturbation. 


Figure  7.2:  Results  obtained  under  acoustic  modulation  in  the  absence  of  flow,  (a)  RMS  pressure  fluctua¬ 
tions  for  acoustic  perturbations  in  phase  (0  =  0).  (b)  RMS  pressure  fluctuations  for  acoustic  perturbations 
in  phase  opposition  (0  =  7 r). 


Case 

N62;i0 

N6f» 

(p=ir 

N6^0 

N67° 

<p=7T 

Ppp/Po 

0 

2.0  % 

0 

2.0  % 

7 r 

1.0  % 

0 

1.0  % 

IT 

Table  7.1:  Test  cases  with  acoustic  modulation. 
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7.2  Calculated  flow  dynamics  and  experimental  data 


When  acoustic  modulations  are  in  phase,  pressure  fluctuations  are  produced  at  the  injector  outlet  at 
the  modulation  frequency.  This  leads  to  oscillations  of  the  axial  velocity  at  the  imposed  frequency. 
Such  variations  in  the  longitudinal  velocity  induce  a  periodic  shedding  of  ring-like  patterns  (Fig.  7.3). 
Structures  are  created  behind  the  inner  lip  and  in  the  outer  mixing  layer.  These  structures  develop  in 
an  axi-symmetric  fashion  producing  annular  vortices.  The  vortices  produced  in  the  inner  jet  lip  entrain 
fluid  from  the  inner  jet,  involving  mass  from  the  inner  jet  in  the  vortex  roll-up.  These  results  are 
qualitatively  consistent  with  experimental  observations  shown  in  Fig.  7.5(a),  where  a  varicose  motion 
is  clearly  identifiable.  When  acoustic  perturbations  are  in  phase-opposition  (Fig.  7.4),  the  resulting 

Hr  jit  -  jfer 

□  □  □ 


Figure  7.3:  Case  “N6^0”  (J=3.05).  Longitudinal  slices  of  density  (white:  60  kg  m  3  ;  black:  maximum 
;  logarithmic  scale).  These  images  correspond  to  three  successive  times  in  one  period. 


transverse  velocity  perturbation  induces  a  sinusoidal  motion,  as  observed  experimentally  in  Fig.  7.5(b) 
and  vortices  are  shed  in  an  asymmetric  fashion. 


Figure  7.4:  Case  (J=3.05).  Longitudinal  slices  of  density  (white:  60  kg  m  3  ;  black:  maximum 

;  logarithmic  scale).  These  images  correspond  to  three  successive  times  in  one  period. 


The  centerline  density  and  axial  velocity  are  plotted  in  Fig.  7.6.  Acoustic  perturbations  significantly 
reduce  the  dense  core  length  and  increase  the  amplitude  of  the  low  velocity  region  just  behind  the  central 
jet.  Characteristic  length  ratios  between  acoustically  perturbed  cases  and  the  reference  case  “N6”  are 
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Figure  7.5:  Experimental  shadowgraphs  for  case  “N6”.  (a)  Acoustic  modulation  is  in  phase,  (b)  Acoustic 
modulation  is  in  phase-opposition. 


given  in  Tab.  7.2.  These  ratios  are  in  agreement  with  experimental  measurements.  It  is  interesting  to 
note  that  the  cases  corresponding  to  phase  opposition  produce  the  shortest  dense  core.  It  is  also  found 
that  multiplication  of  the  perturbation  amplitude  by  a  factor  of  two  does  not  lead  to  a  proportional 
reduction  of  the  dense  core. 


Figure  7.6:  (a)  Density  profile  on  the  centerline,  (b)  Axial  velocity  profile  on  the  centerline. 


Case 

N6™0 

N6^„ 

N6^0 

N6^ 

70.5  //0.5 

1  p  /  Lp,noac 

0.60 

0.34 

0.70 

0.45 

10.1/10.1 

p  / vp,noac 

0.53 

0.36 

0.66 

0.45 

Experiments 

- 

- 

0.68 

0.55 

Table  7.2:  Influence  of  acoustic  modulations  on  the  characteristic  length  in  case  “N6”. 


7.3  Spectral  analysis  of  acoustic  modulation  effects 


The  response  of  the  flow  to  the  external  modulation  may  now  be  discussed  by  examining  the  spectral 
content  of  velocity  fluctuations.  This  is  done  by  plotting  spectral  maps  of  these  fluctuations  along  three 
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axial  lines  in  the  flow.  The  first  (line  C)  corresponds  to  the  centerline,  the  second  (line  B)  is  centered  on 
the  lip,  the  third  (line  A)  is  centered  on  the  outer  jet  annulus  (see  Fig.  7.7). 

Axis  A 
Axis  B 
Axis  C 

□ 

Figure  7.7:  Case  “N6^0”  (J=3.05).  Positions  of  the  three  axial  lines  used  in  the  spatial  spectral  analysis. 


7.3.1  Case  “N6JV 

In  phase  excitation  induces  axial  velocity  modulations  on  the  centerline.  This  is  clearly  seen  in  Fig. 
7.8(a),  7.9(a)  and  7.10(a),  where  the  3  000  Hz  fluctuations  are  well  retrieved.  The  amplitude  of  the 
axial  velocity  fluctuation  is  nevertheless  very  low  in  the  high-density  region  (Fig.  7.11),  in  the  potential 
core  of  the  inner  jet  (x  <1.05  d^,  defined  by  the  region  along  the  centerline  axis  where  p  >  0.99 pi,  see 
Fig.  7.6(a)).  This  region  of  the  flow  is  less  affected  by  the  axial  modulation  than  the  surround  streams. 
The  response  of  the  jets  to  this  modulation  is  determined  by  examining  the  transverse  velocity.  Both 
the  outer  jet  and  the  flow  behind  the  inner  lip  show  a  strong  response  at  3  000  Hz  (Fig.  7.8(b)  and 
7.9(b)).  It  means  that  the  axial  perturbation  promotes  a  destabilization  of  these  jets  and  vortex  roll-up 
at  the  same  frequency.  This  is  not  the  case  along  the  inner  jet  centerline,  where  the  transverse  response 
at  3  000  Hz  is  not  distinguishable  (Fig.  7.10(b)).  Since  the  geometry  and  the  perturbation  are  both 
axisymmetric  around  this  axis,  such  a  velocity  motion  is  naturally  limited.  The  dominant  frequency  after 
the  potential  core  ( x  >1.05  di)  is  close  to  1  000  Hz.  The  spectral  content  of  the  inner  jet  is  shown  in 
Fig.  7.11  at  x  =  di,  at  the  end  of  the  potential  core.  The  axial  velocity  is  dominated  by  the  modulation 
at  3  000  Hz,  with  a  low  amplitude,  but  also  features  a  lower  frequency  component  around  1  200  Hz.  In 
contrast  with  the  outer  jet,  the  radial  velocity  responds  at  half  of  the  imposed  frequency  (around  1  500 
Hz).  These  fluctuations  are  nevertheless  small  compared  to  those  measured  in  the  gaseous  part  of  the 
flow.  As  a  consequence,  the  high  density  part  of  the  inner  jet  seems  to  be  virtually  insensitive  to  the 
acoustic  modulation.  The  frequency  that  dominates  at  the  end  of  the  potential  core  is  close  to  half  of 
the  modulation  frequency.  This  frequency  is  also  observed  in  the  analysis  of  density  fluctuations. 


7.3.2  Case  “N6^w” 

A  transverse  velocity  modulation  is  again  applied  at  3  000  Hz  but  the  phase  is  now  equal  to  7 r.  As  in 
the  previous  case  “N6^0”,  this  transverse  modulation  affects  all  the  jets  as  shown  in  Fig.  7.12  (b),  7.13 
(b)  and  7.14  (b).  The  amplitude  of  this  transverse  perturbation  is  much  lower  in  the  high  density  part  of 
the  inner  jet  (Fig.  7.15  and  7.14(b)  for  x  <  di).  Axial  velocity  perturbations  are  induced  in  the  outer  jet 
and  behind  the  inner  lip  (Fig.  7.12(a)  and  7.13(a)),  associated  with  vortex  roll-up  at  the  same  frequency. 
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Figure  7.8:  Case  “N6^0”  ( J=3.05).  Spectral  map  along  line  A  (r  =  d,/2  +  /,  +  he/2).  (a)  Axial  velocity, 
(b)  Radial  velocity. 
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Figure  7.9:  Case  “N6^0”  (J=  3.05).  Spectral  map  along  line  B  (r  =  di/2-\-li/2).  (a)  Axial  velocity  and 
(b)  Radial  velocity. 
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Figure  7.10:  Case  “N 6^0”  (J=  3.05).  Spectral  map  along  line  C  (r  =  0).  (a)  Axial  velocity,  (b)  Radial 
velocity.  The  dashed  line  represents  the  axial  position  at  which  a  power  spectral  density  is  shown  in  Fig. 
7.11. 
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Figure  7.11:  Case  “N6^0”  (J=  3.05).  Power  spectral  density  of  axial  ( - )  and  radial  ( - )  velocity  in 

the  inner  jet  centerline  at  x  =  d{  (end  of  the  potential  core). 


Again,  the  inner  jet  is  essentially  unaffected,  since  the  measured  amplitudes  are  much  lower  than  for  the 
two  other  cases  (Fig.  7.14(a)).  This  stream  is  dominated  by  an  intermediate  frequency.  A  component  at 
1  500  Hz  is  observed  just  at  the  inner  injector  exit.  This  is  confirmed  on  Fig.  7.15,  where  the  transverse 
modulation  at  3  000  Hz  is  visible,  but  leads  to  a  longitudinal  response  of  the  jet  close  to  1  500  Hz. 
Even  though  this  frequency  is  clearly  identifiable  in  Fig.  7.14(a)  and  7.15,  its  amplitude  is  lower  than 
that  measured  behind  the  lip  or  in  the  outer  jet.  As  previously  observed  for  in  phase  modulation,  the 
sensitivity  of  the  inner  jet  to  the  external  modulation  is  limited.  Again,  the  frequency  which  dominates 
at  the  end  of  the  potential  core  (x  <0.55  d$,  see  Fig.  7.6(a))  is  close  to  half  of  the  modulation  frequency 
and  is  also  observed  in  density  fluctuations. 


5000- 


(a) 


Axis  A 


Axial  velocity 


x/diH 


5000-,, 


(b) 


Radial  velocity 


maximum  =  12  m2  s_1 


Figure  7.12:  Case  (J=  3.05).  Spectral  map  along  line  A  (r  =  di/2  +  U  +  he/ 2).  (a)  Axial 

velocity  and  (b)  Radial  velocity. 


7.4  Discussion 


It  is  found  that  the  inner  jet,  when  submitted  to  an  acoustic  modulation  at  3  000  Hz,  exhibits  fluctuations 
at  the  subharmonic  frequency  of  1  500  Hz.  This  may  be  explained  by  considering  the  stability  properties 
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Figure  7.13:  Case  “NG^r”  {J=  3.05).  Spectral  map  along  line  B  (r  =  di/2  +  k/2).  (a)  Axial  velocity 
and  (b)  Radial  velocity. 
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Figure  7.14:  Case  “NG;^”  (J=  3.05).  Spectral  map  along  line  C  (r  =  0).  (a)  Axial  velocity  and  (b) 
Radial  velocity.  The  dashed  line  represents  the  axial  position  at  which  a  power  spectral  density  is  shown 
in  Fig.  7.15. 


Figure  7.15:  Case  “NG?^”  (J=  3.05).  Power  spectral  density  of  axial  ( - )  and  radial  ( - )  velocity 

on  the  inner  jet  centerline  at  x  =  0.5 di  (end  of  the  potential  core). 
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of  the  central  jet  and  by  analyzing  the  acoustic  transmission  from  the  outer  medium  into  the  high  density 
region. 

It  is  first  interesting  to  situate  the  modulation  frequency  with  respect  to  the  range  of  unstable  fre¬ 
quencies  characterizing  the  jet  and  more  specifically  with  respect  to  the  preferred  frequency  of  the  inner 
jet.  It  is  known  (see  for  example  Birbaud  et  al.  [3]  and  references  in  this  article)  that  when  the  mod¬ 
ulation  frequency  is  well  above  the  preferred  mode  frequency,  the  flow  is  not  affected  by  an  external 
acoustic  modulation.  This  is  best  expressed  in  terms  of  a  Strouhal  numbers  which  can  be  calculated 
by  considering  the  acoustic  modulation  frequency,  the  characteristic  length  and  velocity  associated  with 
each  jet.  The  resulting  “acoustic”  Strouhal  numbers  are  given  in  Tab.  7.3.  While  St^c  and  St^c  are  lower 
or  close  to  the  natural  Strouhal  number  of  the  flow  for  case  “N6”,  St^c  is  higher  than  the  one  measured 
previously  without  modulation.  As  a  consequence,  the  jet  does  not  respond  to  the  modulation  frequency, 
but  oscillates  at  a  subharmonic  frequency  (1  500Hz)  which  is  close  to  the  preferred  frequency  measured 
previously  in  the  absence  of  acoustic  modulation  (Fig.  6.9). 


St  =  ( facl)/u 

length  scale  l 

charac.  vel.  u 

value 

Stic 

he  (0.41  mm) 

ue 

0.087 

Stlac 

li  (0.54  mm) 

ue 

0.115 

stic 

di  (0.51  mm) 

Ui 

0.45 

Table  7.3:  Definition  of  the  Strouhal  numbers  used  in  this  study.  /  represents  the  modulation  frequency 
(3  000  Hz). 


Another  possible  explanation  of  the  absence  of  modulation  at  the  imposed  frequency  may  be  obtained 
by  considering  the  change  in  acoustic  impedance  between  the  outer  stream  and  the  inner  jet.  It  is 
known  that  a  large  impedance  mismatch  reduces  the  transmission  of  waves  across  the  interface.  Acoustic 
impedances  gathered  in  Tab.  7.4  indicate  that  a  large  variation  takes  place  at  the  inner  jet  boundary. 
Under  these  circumstances  acoustic  wave  transmission  into  the  jet  is  reduced  explaining  the  absence  of 
the  modulation  imposed  to  the  system  from  the  outer  region. 


Conditions 

Impedance  (Z  =  pc) 

Chamber 

17  701 

Outer  jet 

19  663 

Inner  jet 

98  700 

Table  7.4:  Specific  acoustic  impedance  in  the  chamber,  inner  and  outer  jets. 


The  previous  arguments  provide  some  clues  on  the  observed  response  of  the  coaxial  flow  to  an  exter¬ 
nally  imposed  transverse  acoustic  motion.  The  strong  reduction  of  the  inner  jet  length  is  certainly  linked 
to  the  strong  response  of  the  outer  jet  to  the  acoustic  modulation,  which  is  not  observed  for  the  inner 
stream.  The  large  vortices  generated  in  the  annular  jet  in  turn  enhance  mixing  between  the  outer  and 
inner  jets,  and  thus  enhances  mixing  of  the  inner  fluid  with  its  surroundings,  leading  to  a  reduction  of 
the  dense  core  length. 
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Chapter  8 

Conclusion 


The  structure  of  non  reactive  transcritical  and  supercritical  jets  is  investigated  by  examining  the  flow 
formed  by  a  coaxial  injector  operating  at  high  pressure.  This  unit  injects  low  temperature  dense  nitrogen 
surrounded  by  an  outer  stream  of  higher  temperature  nitrogen.  The  analysis  relies  on  large  eddy  sim¬ 
ulations  which  are  compared  with  experimental  results.  The  natural  development  of  the  coaxial  flow  is 
examined  for  a  set  of  momentum  flux  ratios.  It  is  shown  that  the  flow  patterns  observed  experimentally 
are  suitably  retrieved  and  that  a  good  correspondence  between  experimental  results  and  simulations  is 
obtained.  It  is  shown  in  particular  that  when  the  momentum  flux  ratio  exceeds  a  critical  value,  the 
central  jet  is  terminated  in  an  abrupt  fashion  at  a  short  distance  from  the  injector  exhaust.  Calculations 
indicate  that  this  is  linked  to  the  formation  of  a  recirculation  region  located  at  a  short  distance  from  the 
injection  section.  This  study  is  complemented  with  a  spectral  analysis  of  the  velocity  components  in  the 
injector  nearfield.  It  is  shown  that  the  frequency  content  is  governed  by  Strouhal  numbers  which  follow 
standard  rules  for  jets  and  wakes.  The  effect  of  an  external  acoustic  modulation  is  then  investigated  by 
imposing  in  phase  and  out  of  phase  perturbations  at  the  boundaries  of  the  domain.  This  is  used  to  study 
the  effect  of  transverse  waves  on  the  dark  core  in  the  jet.  A  strong  reduction  of  the  jet  length  is  clearly 
observed,  both  in  experiment  and  numerics.  The  spectral  content  of  fluctuations  in  the  coaxial  streams  is 
analyzed.  In  the  case  considered  which  corresponds  to  a  Strouhal  number  (based  on  the  inner  jet  velocity 
and  diameter)  about  equal  to  twice  the  Strouhal  number  corresponding  to  the  preferred  mode  of  the  jet, 
it  is  found  that  the  inner  jet  responds  with  a  limited  amplitude  to  the  external  modulation  but  features 
a  subharmonic  motion  at  half  of  the  modulation  frequency.  In  contrast,  the  gaseous  part  of  the  flow  is 
strongly  perturbed  by  the  acoustic  modulation.  Large  scale  vortices  are  released  from  the  outer  annulus 
at  the  modulation  frequency.  It  is  suggested  that  the  relative  insensitivity  of  the  inner  flow  is  linked  to 
the  fact  that  the  frequency  is  above  the  range  of  unstable  frequencies  of  the  central  jet.  It  is  also  found 
that  the  specific  acoustic  impedance  of  the  central  jet  is  much  higher  than  that  corresponding  to  the 
surrounding  stream  and  that  this  may  explain  why  the  modulation  is  not  transmitted  to  the  central  core. 
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N  omenclat  ur  e 


a  Thermal  expansion  coefficient 

/3  Isothermal  compressibility  coefficient 
A  Thermal  conductivity 

p  Dynamic  viscosity 

fil  Turbulent  dynamic  viscosity 

vl  Turbulent  kinematic  viscosity 

uj  Acentric  factor 

p  Density 

p*  Normalized  density  (p  -  Poo)/(pc  -  Poo) 

Poo  Far  field  density 

Poo  Reservoir  initial  density 

pc  Time  averaged  centerline  density 

pe  Outer  jet  injection  density 

pi  Inner  jet  injection  density 

r  Viscous-stress  tensor 

t1  SGS  stress  tensor 

a,  6,  c  Coefficients  for  the  Peng- Robinson  equation  of  state 

cp  Constant  pressure  specific  heat 

cv  Constant  volume  specific  heat 

Cw  WALE  model  constant 

di  Inner  jet  diameter 

E  Total  energy 

es  Internal  energy 
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Spreading  rate  of  the  density  profile 
he  annular  duct  thickness 
hs  Enthalpy 

J  Momentum  flux  ratio 

li  Inner  lip  thickness 

Lp  ’Half  Width  at  Half  Maximum’  of  the  density  profile 


10.1 

V 

^0.1 

XP 

-  x0 

zO.5 

V 

^0.5 

XP 

-  x0 

70.99 

V 

_0.99 

XP 

-  x0 

lu 

Xu  ~ 

lu 

p  Pressure 

Pc  Critical  pressure 

R  Ideal  gas  constant 

r  Specific  gas  constant 

S  Density  ratio 

S  Symmetric  part  of  the  velocity-gradient  tensor 

8  Symmetric  traceless  part  of  the  square  of  the  velocity-gradient  tensor 

T  Temperature 

Tc  Critical  temperature 

Tr  Reduced  temperature 

Reservoir  initial  temperature 
Taver  Duration  of  the  averaging  procedure 
Te  Outer  jet  injection  temperature 
Ti  Inner  jet  injection  temperature 

U  Velocity  ratio 

u  Axial  velocity 

u*  =  u/ue  Normalized  velocity 
Uqq  Reservoir  initial  velocity 

ue  Outer  jet  injection  velocity 

Ui  Inner  jet  injection  velocity 
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x®A  Axial  centerline  position  where  p  =  0.1  {pi  —  pe)  +  pe 

x®‘6  Axial  centerline  position  where  p  =  0.5 (pi  —  pe)  +  pe 

x°p"  Axial  centerline  position  where  p/ pi  =  0.99 

xu  Axial  centerline  position  where  axial  velocity  is  minimum 

q1  SGS  heat  flux  vector 

q  Heat  flux 

v  Velocity  vector 

w  Vector  of  conservative  variables 

Prt  Turbulent  Prandtl  number 

Re  Reynolds  number 

stb  (fue)/(di  +  li) 

Ste  ( fue)/he 

St*  (. fui)/di 

St'  ( fue)/U 
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Appendix  A 

Injection  conditions  for  the 
computed  cases 


Outer  jet  injection  conditions 


Case 

Te  [K] 

me  [g.s-1] 

pe  [kg.m-3] 

ue  [m.s  1] 

Ree 

N2 

152 

1.57 

101 

5.95 

41  600 

N6 

183 

2.69 

73.1 

14.1 

64  800 

N8 

192 

6.46 

68.3 

36.2 

151  000 

Table  A.l:  Injection  conditions  corresponding  to  the  cases  “N2”,  “N6”  and  “N8”  for  the  outer  jet.  m 
denotes  the  mass  flow  rate. 


Inner  jet  injection  conditions 


Case 

Ti  [K] 

rhi  [g.s-1] 

Pi  [kg.m-3] 

ui  [m.s-1] 

Rd 

N2 

117 

0.289 

590 

2.4 

15  000 

N6 

126 

0.292 

420 

3.4 

25  000 

N8 

128 

0.295 

220 

6.6 

52  000 

Table  A. 2:  Injection  conditions  corresponding  to  the  cases  “N2”,  “N6”  and  “N8”  for  the  inner  jet.  m 
denotes  the  mass  flow  rate. 
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Chapter  1 

Abstract 


The  following  work  stems  from  the  spectral  analysis  performed  by  Schmitt  et  al.  [1]  to  the  LES  data 
presented  in  chapters  6  and  7  in  part  I  of  this  report.  In  this  study,  results  from  the  spectral  analysis  of  the 
experimental  data  corresponding  to  two  cases  analyzed  by  Schmitt  et  al  [1]  with  the  LES  techinque  are 
presented  and  compared.  An  introduction  explaining  the  objective  of  spectral  analysis  is  given  in  Chapter 
2  and  the  specific  methods  used  are  provided  in  Chapter  3.  Results  of  this  analysis  are  used  to  provide  a 
better  picture  of  the  behavior  of  transcritical  coaxial  jets  with  and  without  acoustic  excitation.  In  fact, 
spectral  analysis  is  a  technique  that  could  be  used  extensively  on  other  experimental  datasets  already 
obtained  at  AFRL  to  analyze  frequency  content.  In  Chapter  4,  spectral  maps  of  both  experimental  data 
and  simulation  results  are  presented.  It  is  found  that  both  techniques  clearly  identify  the  modulation 
frequency  for  cases  with  acoustic  modulation. 
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Chapter  2 

Introduction 


Spectral  analysis  is  a  useful  tool  when  there  is  a  need  to  obtain  more  information  about  recurring  events 
or  values  contained  in  any  type  of  signal.  From  a  large  continuous  data  set  one  can  extract  the  frequencies 
at  which  these  recurring  events  or  values  take  place.  If  there  is  only  one  recurring  event  on  a  signal,  a 
clear,  strong  mode  will  be  found  using  spectral  analysis.  If  a  signal  has  several  recurring  events  taking 
place  but  there  is  one  that  predominates,  the  frequency  or  mode  of  this  event  will  be  highlighted  over  the 
others  by  this  analysis,  though  it  will  not  appear  in  a  unique  fashion  as  if  there  was  only  one  recurring 
event. 

Spectral  methods  have  been  used  to  determine  whether  a  jet  of  variable  density  issuing  into  ambient 
air  presents  absolute  or  convective  instability  depending  on  the  density  ratio  between  the  jet  and  the 
air  [2].  According  to  the  authors’  analysis,  different  density  ratios  lead  to  different  forms  of  the  power 
spectral  density  (PSD).  One  shows  the  PSD  of  the  jet’s  velocity  with  very  sharp  peaks  in  the  near  field 
of  the  jet  while  the  other  shows  broader  and  less  prominent  peaks.  Using  an  analogy  between  the  PSD 
measured  in  the  wake  of  a  circular  cylinder  (with  prominent  peaks  well  above  the  rest  of  the  signal)  and 
the  PSD  measured  in  the  unstable  shear  layer  of  a  constant  density  isothermal  round  jet  (with  peaks  that 
are  hardly  one  order  of  magnitude  larger  than  the  rest  of  the  signal),  they  conclude  that  the  behavior 
one  of  the  jets  corresponds  to  absolutely  instability  and  the  other  to  convective  instability. 

Raynal  et  al.  [3]  studied  variable-density  jets  using  hot-wire  anemometry.  They  used  their  power 
spectral  results  to  find  that  flow  stability  was  not  affected  by  the  presence  of  the  hot-wire  one  jet  diameter 
downstream  the  exit  of  the  flow  at  several  density  ratios.  They  did  find  that  as  the  probe  was  moved 
downstream  (from  2  to  6  jet  diameters)  the  spectra  was  characterized  by  the  rise  of  a  subharmonic 
frequency  which  is  one  half  the  dominant  frequency.  The  previous  behavior  might  be  associated  with 
the  pairing  process  of  periodic  structures  present  in  the  flow.  They  also  observed  an  increase  of  the 
spectral  level,  which  reflects  the  spatial  amplification  of  the  perturbations  as  the  hot-wire  was  moved 
downstream.  The  power  spectra  also  showed  that  when  the  jet- to- ambient  density  ratio  is  increased  from 
0.31  to  0.61  the  peak  frequency  decreases  several  orders  of  magnitude  until  it  reaches  the  same  size  of 
the  other  smaller  peaks  at  0.51  and  reaches  the  same  magnitude  of  the  background  noise  at  0.61. 

In  their  work  identifying  combustion  instabilities  in  a  multiple  inlet  combustor,  Poinsot  et  al.  [4] 
use  the  power  spectral  density  maps  of  the  signal  detected  by  a  microphone  and  plot  it  as  a  function 
of  air  flow  rate  and  equivalence  ratio.  Each  map  spans  a  frequency  band  centered  around  one  resonant 
frequency  of  the  system.  They  found  that  each  mode  occurs  in  a  restricted  domain  of  parameters  and  as 
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these  parameters  are  increased  the  trend  is  to  have  higher  acoustic  intensities  and  frequencies. 

For  this  brief  analysis,  we  will  use  spectral  analysis  as  a  tool  to  identify  periodicities  in  our  data  sets. 
These  sets  have  been  obtained  via  experiments  and  computer  simulations  of  a  cryogenic,  transcritical 
coaxial  jet.  On  the  experimental  side,  the  signal  is  the  intensity  of  light  at  each  pixel  from  the  images 
captured  by  a  high-speed  camera.  On  the  modeling  side  our  data  consists  of  the  values  of  the  different 
variables  at  a  given  location  or  node  in  the  simulation  mesh.  There  are  two  cases  which  are  studied  both 
experimentally  and  computationally.  The  effects  of  acoustics  are  also  studied  for  one  of  the  cases.  These 
two  cases  will  be  compared  and  the  results  will  be  shown. 
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Chapter  3 

Methods,  assumptions  and 
procedures 


The  experimental  setup  and  the  modeling  method  are  the  same  as  those  described  in  sections  3  and  5  in 
part  I  of  this  report.  For  reference,  the  injector  used  in  the  experimental  study  and  the  computational 
model  with  respective  lengths  are  reproduced  in  Figs.  3.1  and  3.2  below. 


Figure  3.1:  (a)  Photograph  of  the  shear  coaxial  injector  used  in  this  study,  (b)  Exit  plane  view  of  the 
injector  with  characteristic  dimensions. 


In  order  to  analyze  the  spectral  content  of  the  data  Welch’s  method  was  used  in  this  study.  This 
method  allows  the  user  to  obtain  an  estimate  of  the  intensity  of  the  signal  as  a  function  of  frequency. 
Welch’s  method  is  based  on  the  method  of  averaged  periodograms  also  known  as  Bartlett’s  method  [5]. 
Periodograms  are  used  to  estimate  the  spectral  content  of  a  signal,  but  they  present  inherent  problems 
such  as  large  variance  of  the  data  and  measurement  noise.  Bartlett’s  method  takes  the  signal  and  divides 
it  in  smaller  segments,  then  computes  the  periodogram  of  each  segment  and  finally  finds  the  average  of 
the  periodograms.  This  helps  decrease  variance  issues  with  the  data  in  exchange  for  the  ability  to  resolve 
smaller  frequency  ranges  [6]. 
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Welch’s  method  also  uses  the  same  type  of  averaging  introduced  by  Bartlett’s  method  and  adds 
a  procedure  that  helps  reduce  the  noise.  First  the  signal  is  divided  in  overlapping  segments  and  the 
amount  of  overlapping  might  be  changed  as  desired  by  the  user.  If  the  all  the  data  points  in  the  signal 
are  part  of  two  different  segments  the  overlapping  is  maximized.  This  is  referred  to  as  a  50  %  overlap. 
If  there  is  no  overlap  assigned  by  the  user  one  returns  to  Bartlett’s  method.  After  creating  overlapping 
segments,  the  segments  will  be  altered  by  a  so-called  window  function  which  helps  reduce  the  noise  but 
generally  gives  more  weight  to  the  data  at  the  center  of  the  segment  and  less  to  the  data  at  the  edges, 
thus  the  reason  for  the  data  overlap.  Once  the  segments  are  overlapped  and  “windowed”  the  periodogram 
of  each  is  found  and  then  the  periodograms  are  averaged.  The  results  is  the  graph  of  intensity  or  power 
spectral  density  of  the  signal  as  a  function  of  frecuency  [6,  7,  8]. 


Case  2  (  J  ~  1.05,  Pr  =  1.05  ) 

T  [K] 

m  [g.s  ■*■] 

p  [kg.m  3] 

u  [m.s  ■*■] 

Re 

outer  jet 

152 

1.57 

101 

5.95 

41  600 

inner  jet 

117 

0.289 

590 

2.4 

15  000 

Table  3.1:  Injection  conditions  for  the  outer  jet  and  for  the  inner  jet  corresponding  to  the  cases  “N2” 
(simulation)  &  “near2”  (experiment).  J  is  the  momentum  flux  ratio,  Pr  is  the  reduced  pressure  and  m 
denotes  the  mass  flow  rate. 


Case  6  (  J  ~  3.0,  Pr  =  1.05  ) 

T[  K] 

rh  [g.s  -1] 

p  [kg.m  3] 

u  [m.s  1] 

Re 

outer  jet 

183 

2.69 

73.1 

14.1 

64  800 

inner  jet 

126 

0.292 

420 

3.4 

25  000 

Table  3.2:  Injection  conditions  for  the  outer  jet  and  for  the  inner  jet  corresponding  to  the  cases  “N6” 
(simulation)  &  “near6”  (experiment).  J  is  the  momentum  flux  ratio,  Pr  is  the  reduced  pressure  and  m 
denotes  the  mass  flow  rate. 

For  the  simulation  data  the  signals  analyzed  were  radial  velocity,  axial  velocity  and  pressure.  Spectral 
graphs  representing  the  power  spectral  density  of  these  signals  as  a  function  of  frequency  at  a  given 
downstream  location  from  the  inner  jet  exit  were  obtained.  Spectral  maps  were  built  showing  contours 
of  the  signal  power  spectral  density  at  different  frequencies  for  all  downstream  locations.  The  graphs  and 
maps  obtained  from  simulation  data  are  shown  in  Sections  6.3  and  7.3  in  part  I  of  this  report  and  most 
of  them  are  reproduced  in  this  part  of  the  report  for  comparison. 
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In  the  case  of  experimental  measurements,  the  signal  consisted  of  the  intensity  of  light  in  the  recorded 
images.  The  same  procedure  followed  to  extract  the  power  spectral  density  from  simulation  data  was 
used  with  the  experimental  results.  Spectral  maps  were  also  built  showing  contours  of  the  signal  power 
spectral  density  at  different  frequencies  for  all  downstream  locations.  In  order  to  compare  results  with 
experiments,  the  same  axes  selected  to  analyze  the  simulation  data  were  selected  for  the  spectral  analysis 
of  the  experimental  data.  An  image  showing  the  axes  superimposed  to  a  shadowgraph  of  the  coaxial  jet 
is  shown  in  Fig.  3. 3. a  while  the  image  shown  to  illustrate  the  axes  in  section  7.3  in  part  I  of  this  report 
is  reproduced  in  Fig.  3.3.b  below. 


Axis  A 
Axis  B 
Axis  C 


Figure  3.3:  Positions  of  the  three  axial  lines  used  in  the  spatial  spectral  analysis,  (a)  Experimental  case 
“near6”.  (b)  Simulation  case  “N6^0”. 


To  apply  Welch’s  method  both  the  simulation  and  experimental  data  sets  were  segmented  in  M 
overlapping  blocks  with  a  50  %  overlap  to  minimize  variance  problems.  Because  the  number  of  samples 
is  limited  the  averaging  is  carried  out  over  M  =  8  blocks  in  the  outer  flow  region  and  on  M  =  4  blocks  in 
the  inner  jet  region  for  the  simulations  and  over  M  =  8  blocks  for  all  experimental  data.  The  frequency 
resolution  corresponding  to  the  different  locations  along  the  axes  shown  in  Fig.  3.3  is  given  in  Table  3.3. 


St  =  (. fl)/u 

l  (Fig.  3.1.b) 

u 

Axis  (Fig.  3.3) 

A  /  [Hz] 

St 1 

di  (0.51  mm) 

Ui 

c 

N2:  187;  N6:  125;  near2  &  near6:  200 

St1 

li  (0.54  mm) 

lie 

B 

N2:  375;  N6:  250;  near2  &  near6:  200 

Ste 

he  (0.41  mm) 

Ue 

A 

N2:  375;  N6:  250;  near2  &  near6:  200 

Table  3.3:  Definition  of  characteristic  Strouhal  numbers.  The  characteristic  length  and  velocity  l  and  u 
are  given  in  the  second  and  third  columns.  The  frequency  /  corresponds  to  the  power  spectral  density 
maximum  from  simulations  and  experiments.  A /  is  the  frequency  resolution. 
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Chapter  4 

Results  and  discussion 

4.1  Spectral  analysis  of  cases  without  acoustic  modulation 


For  cases  without  acoustic  modulation,  it  is  of  interest  to  identify  any  periodicities  in  the  flow  dynamics. 
The  spectral  maps  of  these  cases  thus  shed  some  light  into  the  natural  modes  that  might  be  present  in 
the  coaxial  jet’s  dynamics.  This  information  is  useful  when  identifying  modes  that  might  be  enhance  by 
the  modulation  frequency  or  to  find  approximate  Strouhal  numbers  for  the  natural  modes  present.  In 
the  following  figures  the  spectral  maps  of  transcritical  cases  “near2”  and  “near6”  are  shown.  The  plots 
show  frequency  in  Hertz  versus  distance  downstream  the  coaxial  injector’s  exit  in  inner  jet  diameters  (di) 
where  one  di  corresponds  to  0.51  millimeters. 


Figure  4.1:  Spectral  maps  of  experimental  cases  “near2”  and  “near6”  without  acoustic  modulation  along 
the  inner  jet  center  axis  or  axis  C  as  shown  in  Fig.  3.3.  (a)  Case  unear2”.  (b)  Case  “near6”. 


In  the  spectral  maps  we  can  observe  the  different  values  that  the  intensity  of  the  signal  reaches.  In  Fig. 
4. 5. a  along  the  central  axis  we  observe  a  mode  of  approximately  1.2  kHz  towards  the  end  of  the  image 
(around  25  di  downstream)  for  case  near2,  this  might  appear  as  a  strong  mode  but  its  intensity  is  low 
compared  to  the  maximum  value  in  Fig.  4.5.b.  This  map  presents  modes  near  1.3  and  2.1  kHz  between 
8  to  10  di  downstream.  The  location  of  the  modes  follows  the  same  trend  as  the  length  of  the  inner  jet 
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core  with  a  lower  momentum  flux  ratio  case  such  as  near2  having  a  mode  present  further  downstream  as 
compared  to  the  higher  momentum  flux  ratio  case  near6 


Figure  4.2:  Spectral  maps  of  experimental  cases  “near2”  and  “near6”  without  acoustic  modulation  along 
the  center  of  the  inner  jet  post  or  axis  B  as  shown  in  Fig.  3.3.  (a)  Case  “near2”.  (b)  Case  “near6”. 


The  maps  in  Fig.  4.8  show  different  dynamics  then  the  ones  in  the  previous  image.  Now  the  maximum 
intensity  of  the  signal  for  both  cases  is  similar  and  the  low  momentum  flux  ratio  case  shows  a  mode  near 
1.7  kHz  between  6  and  8  di  downstream  whereas  the  higher  momentum  flux  ratio  case  shows  a  mode 
near  2.8  kHz  also  around  7  and  8  di  downstream.  A  natural  mode  near  3  kHz  is  always  interesting  since 
modulation  has  been  applied  at  this  frequency  and  certain  phenomena  could  be  better  understood  if  we 
find  there  is  coupling  between  a  natural  mode  and  the  imposed  acoustic  forcing. 


Figure  4.3:  Spectral  maps  of  experimental  cases  “near2”  and  “near6”  without  acoustic  modulation  along 
the  outer  jet  center  axis  or  axis  A.  as  shown  in  Fig.  3.3.  (a)  Case  “near2”.  (b)  Case  “near6”. 


For  the  last  two  maps  in  Fig.  4.5  we  focus  on  the  dynamics  of  the  outer  jet  axis.  We  notice  that  the 
maximum  intensity  of  the  signal  for  case  “near2”  is  higher  than  any  of  the  other  maps  we  have  seen  so 
far.  The  mode  with  such  intensity  is  located  15  di  downstream  the  inner  jet  exit  and  at  a  frequency  of 
1.2  kHz.  More  interestingly,  the  case  “near 6”  features  a  region  of  maximum  intensity  very  close  to  3.0 
kHz;  however,  the  rest  of  the  map  shows  regions  of  high  intensity  at  several  locations  and  frequencies, 
which  make  a  clear  mode  difficult  to  distinguish. 
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Figure  4.4:  Power  spectral  density  of  radial  velocity  for  simulation  data,  (a)  Case  “N2”  in  the  mixing 
layer  behind  the  inner  lip  (this  spectrum  is  found  3  di  downstream  along  the  B  axis  as  shown  in  Fig. 
3.3).  (b)  Case  “N6”  in  the  inner  jet  centerline  (this  spectrum  is  found  2  di  downstream  along  the  C  axis 
as  shown  in  Fig.  3.3). 


The  two  graphs  shown  in  Fig.  4.4  show  results  without  acoustic  forcing  found  from  simulation  data 
that  have  similar  trends  to  those  obtained  from  the  experimental  images.  The  first  graph  (Fig.  4. 4. a) 
shows  a  distinctive  wide  peak  at  a  location  along  the  center  of  the  inner  jet  post  (B  axis)  near  1.6  kHz  for 
the  simulation  case  “N2”  which  is  the  equivalent  to  the  experimental  case  “near 2”  which  shows  a  mode 
which  is  not  as  wide  but  still  strong  near  this  frequency  too  in  Fig.  4. 8. a.  The  main  difference  however 
is  the  physical  location  of  the  mode.  In  the  simulation  we  find  the  mode  at  a  location  3  di  downstream 
while  the  experimental  results  show  a  mode  predominantly  between  6  to  10  di  downstream  the  injector 
exit. 

In  Fig.  4.4.b  there  are  two  peaks  that  are  clearly  visible  near  0.95  and  1.7  kHz.  These  peaks  are  close 
to  the  same  modes  found  in  the  spectral  map  of  Fig.  4.5  where  we  observe  maximum  intensities  near 
1.3  and  2.1  kHz.  Again,  for  the  simulation  we  find  the  location  of  these  nodes  at  2  di  downstream  while 
the  experimental  ones  are  located  between  8  to  10  di  downstream.  Since  the  location  of  the  probe  for 
the  simulation  case  is  related  to  the  end  of  the  potential  core  as  explained  in  section  6.3.1  in  part  I,  the 
phenomena  in  both  cases  might  be  equivalent  since  the  end  of  the  dark  core  for  the  experimental  case 
was  measured  at  6  di  downstream  the  injector. 
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4.2  Spectral  analysis  of  a  case  with  acoustic  modulation 


Spectral  maps  of  experimental  case  “near6”  and  simulation  case  “N6”  with  acoustic  modulation  are 
presented  in  this  section.  The  spectral  information  shows  a  clear  influence  of  the  acoustic  modulation  with 
maps  with  very  large  maximum  intensities  and  modes  that  correspond  to  the  excitation  frequency.  This 
is  further  proof  of  the  strong  effect  that  the  acoustic  modulation  has  on  the  flow  even  if  the  perturbation 
is  only  a  small  fraction  (1%  or  2%  maximum)  of  the  mean  values. 


Figure  4.5:  Spectral  maps  of  experimental  case  “near6”  with  acoustic  modulation  along  the  inner  jet 
centerline  or  axis  C  as  shown  in  Fig.  3.3.  (a)  In-phase  modulation,  (b)  Out-of-phase  modulation. 
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Figure  4.6:  Spectral  maps  of  axial  velocity  along  the  inner  jet  centerline  or  axis  C  (r  =  0)  for  (a)  Case 
“N6^0”  (in-phase  modulation)  and  (b)  Case  (out-of-phase  modulation). 


In  the  first  set  of  maps  the  dynamics  along  the  central  axis  are  shown.  For  the  case  in  Fig.  4. 5. a 
with  acoustic  modulation  in  phase  a  mode  just  over  3.0  kHz  can  be  observed  between  6  and  7  di.  This  is 
basically  the  same  mode  observed  for  the  case  with  out-of-phase  modulation  in  Fig.  4.5.b.  Interestingly, 
the  modes  are  located  at  the  location  downstream  where  the  inner  jet  dark  core  ends  when  no  acoustic 
modulation  is  applied.  This  means  that  the  zone  where  the  last  dense  packets  of  fluid  exist  when  no 
acoustic  excitation  is  applied  now  becomes  a  region  where  intense  flow  periodicities  and  mixing  take 
place.  In  fact,  if  one  looks  at  the  spectral  maps  from  the  simulation  data,  on  observes  modes  at  the 
forcing  frequency  of  3  kHz  that  reach  a  maximum  5  di  downstream  for  the  axial  velocity  signal  when  the 
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Figure  4.7:  Spectral  maps  of  radial  velocity  along  the  inner  jet  centerline  or  axis  C  (r  =  0)  for  (a)  Case 
“N6^0”  (in-phase  modulation)  and  (b)  Case  “NG^L^”  (out-of-phase  modulation). 


modulation  is  in  phase  (Fig.  4. 6. a)  and  for  the  radial  velocity  signal  when  modulation  is  out  of  phase 
(Fig.  4.7.b).  Both  images  clearly  show  a  mode  that  might  extend  downstream  beyond  the  5  di  shown 
which  would  correspond  very  closely  with  the  mode  observed  in  Fig.  4.5  for  the  experimental  cases. 


Figure  4.8:  Spectral  maps  of  experimental  case  “near6”  with  acoustic  modulation  along  the  center  of  the 
inner  jet  post  or  axis  B  as  shown  in  Fig.  3.3.  (a)  Modulation  in  phase,  (b)  Modulation  out  of  phase. 


The  spectral  map  in  Fig.  4. 8. a  also  shows  a  mode  at  the  same  frequency,  slightly  over  3.0  kHz,  as 
the  one  shown  in  the  previous  images.  This  is  expected  since  the  modulation  applied  was  the  same  for 
all  conditions.  The  interesting  feature  of  this  image  is  the  “distance  of  influence”  of  the  mode  which 
extends  from  3  to  7  di  downstream.  Such  a  region  could  indicate  that  a  larger  “periodic  mass  transfer” 
or  recirculation  zone  is  created  just  after  the  inner  jet  post.  For  the  out-of-phase  condition  in  Fig. 
4.8.b,  however,  no  such  band  exists  but  instead  there  is  a  point,  which  highlights  a  mode  with  a  level  of 
intensity  an  order  of  magnitude  larger  than  its  in-phase  counterpart.  This  very  strong  mode  shares  the 
same  frequency  though  it  is  located  between  2  and  3  di  downstream.  This  could  be  the  consequence  of 
strong  transversal  velocity  oscillations  acting  near  the  exit  of  the  injector  and  quickly  mixing  the  dense 
inner  jet  shortly  afterwards.  An  interesting  feature  of  these  last  two  spectral  maps  has  to  do  with  the 
region  where  the  strongest  modes  appear.  Up  to  this  point  no  previous  experimental  spectral  map  had 
shown  any  significant  activity  in  the  region  between  0  and  5  di  downstream  except  maybe  for  the  case 
“near  6”  along  the  outer  jet  centerline  in  Fig.  4.11.b  where  there  is  some  activity  in  the  sub  1kHz  range 
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Figure  4.9:  Spectral  maps  of  axial  velocity  along  the  inner  jet  post  centerline  or  axis  B  (r  =  di/2  +  k/ 2) 
for  (a)  Case  “N6^0”  (in-phase  modulation)  and  (b)  Case  “N6|'y ”  (out-of-phase  modulation). 
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Figure  4.10:  Spectral  maps  of  radial  velocity  along  the  inner  jet  post  centerline  or  axis  B  (r  =  g^/2  +  ^/2) 
for  (a)  Case  “N6^0”  (in-phase  modulation)  and  (b)  Case  (out-of-phase  modulation). 
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and  more  recently  in  our  out-of-phase  modulation  case  “near2”  along  the  inner  jet  centerline  in  Fig.  4.5.b 
where  the  3kHz  mode  starts  showing  around  4  di  downstream.  Now  we  have  two  very  strong  modes, 
specially  the  one  present  when  modulation  is  out  of  phase  in  the  region  between  the  two  jets  shortly 
before  the  jets  exit  the  injector  and  start  mixing.  To  confirm  this  new  trend,  the  spectral  maps  can  be 
compared  with  the  spectral  data  from  simulations.  The  similarities  could  not  be  closer  since  the  radial 
velocity  signal  for  the  in-phase  modulation  case  (4. 10. a)  shows  a  relatively  broad  mode  starting  just  after 
1  di  and  continuing  until  3  di  and  extending  weakly  into  5  di  downstream;  however,  when  the  radial 
velocity  mode  loses  strength  the  axial  velocity  one  predominates  as  seen  in  Fig.  4. 9. a  and  extends  until 
the  end  of  the  map.  This  “combined”  broad  mode  is  similar  to  the  broad  mode  observed  in  experimental 
case  in  Fig.  4. 8. a.  The  same  comparison  between  simulation  and  experiment  could  be  drawn  by  observing 
the  case  “N6”  under  out-of-phase  modulation,  where  the  combined  regions  of  influence  of  the  axial  and 
radial  velocity  modes  in  Figs.  4.10.b  and  4.10.b,  result  in  a  larger  mode  extending  from  approximately 
0.5  to  2.5  di  downstream.  This  is  a  very  close  match  to  the  region  of  influence  of  the  mode  found  from 
case  “near6”  under  out-of-phase  modulation  and  shown  in  Fig.  4.8.b. 


Figure  4.11:  Spectral  maps  of  experimental  case  “near6”  with  acoustic  modulation  along  the  outer  jet 
center  axis  or  axis  A.  as  shown  in  Fig.  3.3.  (a)  In-phase  modulation,  (b)  Out-of-phase  modulation. 
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Figure  4.12:  Spectral  maps  of  axial  velocity  along  the  outer  jet  flow  centerline  or  axis  A  (r  =  di/2  +  k  + 
he/2)  for  (a)  Case  “N6^0”  (in-phase  modulation)  and  (b)  Case  (out-of-phase  modulation). 


The  last  two  spectral  maps  are  shown  in  Fig.  4.11  and  they  represent  the  center  region  of  the  annular 
jet.  Both  conditions  present  phenomena  at  the  excitation  frequency  of  3  kHz.  The  condition  with  in-phase 
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Figure  4.13:  Spectral  maps  of  radial  velocity  along  the  outer  jet  flow  centerline  or  axis  A  (r  =  di/2  +  //  + 
he/ 2)  for  (a)  Case  “N6?=0”  (in-phase  modulation)  and  (b)  Case  (out-of-phase  modulation). 


modulation  shows  again  a  mode  that  is  spread  over  a  few  di  but  the  mode  is  farther  downstream  than  the 
one  in  Fig.  4. 8. a.  The  mode  is  thus  now  located  between  8  and  11  di  and  with  similar  overall  strength 
but  now  with  50%  larger  maximum  intensity.  Similarly,  the  out-of-phase  condition  presents  a  very  well 
defined  mode  as  that  shown  in  Fig.  4.8.b  but  now  between  6  and  7  di  downstream.  This  trend  might 
give  us  an  idea  of  the  direction  the  strongest  interactions  follow  as  they  are  convected  downstream.  Also 
wider  spread  of  the  modes  along  the  interrogation  line  is  observed  when  the  modulations  are  in-phase 
for  axes  A  and  B  as  compared  to  the  spread  of  the  modes  when  the  modulations  are  out-of-phase.  The 
periodicities  in  the  outer  layers  of  the  flow  are  apparently  more  localized  when  the  flow  has  large  velocity 
fluctuations  at  its  center.  The  simulation  results  in  Fig.  4.12  support  the  already  expected  findings  that 
strong  modes  observed  at  the  modulation  frequency  and  they  also  show  combined  broad  modes;  however, 
their  locations  differ,  since  the  nodes  from  the  simulation  spectral  maps  show  stronger  nodes  between 
2  and  4  di  downstream  and  the  spectral  maps  from  experiments  show  these  modes  in  the  6  to  10  di 
downstream  zone. 
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Chapter  5 

Conclusions 


The  study  of  the  spectral  content  of  the  experimental  images  is  useful  to  further  support  the  findings 
drawn  from  direct  image  processing  and  simulation  results.  This  analysis  provides  another  tool  to  inter¬ 
rogate  the  behavior  of  the  flow.  It  is  an  effective  manner  to  identify  periodicities  which  otherwise  could 
go  unnoticed.  For  instance,  the  low  momentum  flux  ratio  case  studied,  “near2”  featured  modes  in  the  1 
to  2  kHz  range  at  different  locations  downstream  the  jet  exit  when  no  acoustic  modulation  was  applied. 
This  is  an  interesting  information  when  trying  to  decide  if  the  imposed  frequency  during  acoustic  forcing 
might  be  exciting  existing  modes.  This  could  make  a  difference  in  previous  experimental  findings  showing 
that  small  acoustic  perturbations  lead  to  large  dark  core  length  reductions  and  enhanced  mixing,  since  it 
would  not  be  the  same  it  the  perturbations  are  imposed  at  a  resonant  frequency  or  mode  of  the  coaxial 
jet. 


Another  instructive  feature  of  this  analysis  is  the  close  agreement  between  experimental  and  simulation 
data  when  comparing  spectral  maps  of  the  cases  “near6”  and  “N6”  with  acoustic  modulation  even  when 
the  signals  compared  were  not  the  same  (light  intensity  from  each  pixel  for  the  experimental  data  and 
radial  and  axial  velocity  for  the  simulation  data).  The  expected  forcing  frequency  of  3kHz  was  clearly 
present  in  all  spectral  maps  obtained  from  experimental  data  and  these  findings  were  in  agreement  with 
the  spectral  maps  from  the  simulation  data,  not  only  qualitatively  but  also  qualitatively  in  regard  to  the 
distance  downstream  from  the  injector  at  which  the  modes  were  found.  These  is  more  significant  given 
the  fact  that  the  spectral  information  of  the  non- acoustic  case  did  not  agree  to  the  same  extent.  On 
the  other  hand,  the  comparison  between  the  experimental  data  and  simulation  data  for  the  case  without 
acoustics  was  not  as  extensive. 

Further  analysis  of  experimental  and  simulation  data  could  be  performed  to  identify  periodicities 
and  relate  them  to  the  observed  features  of  the  flow.  The  spectral  data  could  be  useful  to  identify 
characteristic  modes  of  the  coaxial  jet  and  obtain  the  associated  non-dimensional  frequencies  (Strouhal 
numbers).  Then  a  map  showing  the  trends  of  these  modes  as  various  coaxial  jet  parameters  are  changed 
could  be  created.  This  could  help  elucidate  certain  types  of  behavior  observed  in  recent  experimental 
results,  such  as  the  inversion  of  the  effects  of  in-phase  versus  out-of-phase  modulation  for  different  types 
of  injectors  [9]  and  extreme  reduction  and  even  obliteration  of  the  inner  jet  core  when  forced  at  certain 
frequencies  [9,  10]. 
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N  omenclat  ur  e 


di  inner  jet  diameter 

A /  resolution  frequency  of  spectral  data 
e  subscript  or  superscript  for  outer  jet  flow 

he  outer  jet  thickness 

i  subscript  or  superscript  for  inner  jet  flow 

/  frequency,  characteristic  frequency 

Hz  Hertz  or  cycles  per  second  [m-1] 

J  Momentum  flux  ratio  (J  =  peu^/piuf) 

kHz  kiloHertz  or  1000  cycles  per  second  [103  m-1] 

l  characteristic  length,  subscript  for  inner  jet  post  or  lip 

LES  Large  Eddy  Simulation 

li  inner  jet  post  thickness 

M  Number  of  blocks  used  for  spatial  averaging 

m  mass  flow  rate 

simulation  case  N6  with  2%  amplitude  and  in-phase  modulation 

simulation  case  N6  with  2%  amplitude  and  out-of-phase  modulation 
P  pressure 

Pc  critical  pressure 

Pr  reduced  pressure  (Pr  =  P/Pc ) 

PSD  Power  Spectral  Density 

r  radial  (transversal)  distance  from  the  center  axis  of  the  flow 

Re  Reynolds  Number 

p  density 

St  Strouhal  Number 

T  Temperature 

u  velocity,  characteristic  velocity 
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Direct  Numerical  Simulation  of 
Flows  under  Transcritical  Conditions 
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Chapter  1 

Abstract 


The  purpose  of  this  study  is  to  develop  and  validate  a  direct  numerical  simulation  (DNS)  code  used  to 
analyze  the  behavior  of  transcritical  flows.  The  code,  which  is  intended  to  become  a  general  purpose  DNS 
code  with  an  emphasis  on  combustion,  has  been  developed  by  the  numerical  combustion  group  at  the 
EM2C  laboratory  of  Ecole  Centrale  Paris.  It  comprises  different  thermodynamic  modules  which  range 
from  simple  air  calculations  without  combustion  to  complex  real  gas  thermodynamics  of  flows  close  to 
the  critical  point  of  the  substance.  An  initial  validation  of  this  code  is  first  carried  out  by  observing  the 
evolution  of  an  initial  temperature  profile  with  a  basic  Gaussian  shape.  This  validation  process  provides 
the  order  of  the  numerical  scheme.  The  code  is  next  used  to  test  situations  that  will  be  of  interest  to  the 
transcritical  evolution  of  shear  flows  in  parts  I  and  II  of  this  report.  The  case  of  a  two-dimensional  mixing 
layer  is  envisaged.  Initial  calculations  are  carried  out  in  the  case  of  a  two-dimensional  (2D)  problem  in 
which  two  streams  form  an  initial  vortex  sheet.  The  results  of  this  analysis  are  used  to  provide  a  picture 
of  how  DNS  tools  can  be  used  to  study  transcritical  behavior  in  these  types  of  flows. 
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Chapter  2 

Introduction 


Direct  Numerical  Simulation  (DNS)  is  a  useful  modeling  tool  when  the  geometry  of  the  problem  is  not 
complex  and  the  difference  between  the  largest  and  smallest  scales  of  the  flow  is  not  too  large.  In  practice, 
very  few  applications  have  the  characteristics  mentioned  above  so  DNS  is  often  used  to  give  the  researcher 
the  opportunity  to  study  more  in  depth  one  characteristic  of  a  more  complex  problem  or  to  help  develop 
models  that  can  be  later  used  in  Large  Eddy  Simulations  (LES)  to  handle  complex  geometries  such  as 
those  required  in  aircraft  and  space  propulsion  problems. 

Research  involving  DNS  for  supercritical  flows  has  been  done  by  Okong’o  and  Bellan  [1,  2]  on  the 
characteristics  of  high  pressure  transitional  flows.  They  formulate  their  equations  according  to  fluctuation- 
dissipation  theory  which  adds  the  thermal  diffusion  factor  to  the  list  of  transport  properties  playing  a 
role  in  high  pressure  systems.  They  observe  that  dissipation  is  mostly  due  to  species  mass  flux  and 
turbulence  models  should  focus  in  this  property  rather  than  the  momentum  flux.  Selle  et  al  [4]  analyze 
this  transitional  supercritical  mixing  layer  in  detail  with  DNS  to  observe  the  behavior  of  the  small  scales 
of  the  flow  and  compare  that  to  the  subgrid-scale  (SGS)  closure  models  used  in  LNS  techniques.  Selle 
and  Ribert  also  investigate  the  impact  of  the  equation  of  state  used  in  DNS  on  turbulence  models  for 
supercritical  flows  [3]. 

At  the  EM2C  laboratory  at  Ecole  Centrale  Paris,  the  Numerical  Combustion  Group  has  developed  a 
code  nicknamed  YWC  (Yes  We  Can)  which  uses  DNS  to  investigate  different  types  of  problems.  The  code 
is  fully  parallel  and  has  many  options  and  functionalities  that  can  be  set  up  by  the  user.  Among  them  are 
the  posibility  of  performing  a  calculation  in  one,  two  or  three  dimensions.  The  type  of  numerical  scheme 
used  (centered,  compact)  and  the  possibility  to  use  filtering.  The  user  also  can  select  which  boundary 
conditions  to  use:  periodic,  inlet,  outlet,  wall,  adiabatic  wall  and  even  there  is  an  option  to  create  a 
special  boundary  condition  if  a  more  complex  setup  is  needed. 

In  particular,  the  code  YWC  features  different  thermodynamic  modules  to  solve  various  types  of 
problems.  The  simplest  is  the  module  AIR,  which  only  deals  with  this  fluid  (air)  as  a  single  species  and 
for  cases  with  no  chemical  reactions.  A  second  module,  SIMPLE  CHEMISTRY,  deals  with  reactions 
considering  only  a  few  species.  This  simple  chemistry  option  allows  for  faster  calculation  of  problems 
where  the  details  of  the  reaction  are  not  an  issue.  The  other  three  modules  deal  with  either  perfect  or  real 
gas  thermodynamic  problems.  The  basic  one,  REGATH,  features  more  species  and  handles  more  complex 
combustion  problems.  A  “light”  version  called  REGTAH  LITE  was  developed  to  decrease  calculation 
costs  when  not  all  the  functionality  of  REGATH  was  needed.  A  final  module,  REGATH  TRC  was 
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conceived  for  problems  dealing  with  fluids  near  their  critical  point.  It  is  this  last  module  which  is  better 
suited  to  study  transcritical  flow  problems  such  as  the  ones  described  in  this  report. 
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Chapter  3 

Methods,  assumptions  and 
procedures 

3.1  YWC  governing  equations 


The  direct  numerical  simulation  tool  YWC  solves  the  balance  of  mass,  species,  momentum  and  energy 
equations.  These  are  partial  differential  equations  which  require  extensive  computational  capability  if 
direct  calculation  is  desired,  especially  for  flows  with  large  differences  between  their  largest  and  smallest 
scales  (high  Re  flows). 


Conservation  of  mass 

dp  dpUi 

dt  dxi 


Conservation  of  species 

dpYk 

dt 


dpUiYk  d 

- b - 1"  Uk  +  ~r\ — 

dxi  dxi 


where  subscript  k  refers  to  the  species  (k  G  [1,7V]) 


Conservation  of  momentum 

dpUj  dpUiUj  dp  di~ij 

dt  dxi  dxj  +  dxi  ^  ^ 

Tij  is  the  (i,j)  element  of  the  viscous  tensor  r,  which  is  calculated  by  : 


Tij  =  p 
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Conservation  of  energy 


dpe  dpuie  dpui  dr^Uj  dqi 

dt  dxi  dxi  dxi  dxi 


(3.14) 


3.2  YWC  numerical  aspects 

3.2.1  Temporal  discretization 


A  numerical  solver  needs  to  discretize  a  domain  in  space  and  also  in  time.  The  method  used  to  discretize 
the  time  domain  in  the  YWC  solver  is  the  4th  order  Runge-Kutta  scheme.  For  the  following  differential 
equation  : 


dp 

dt 


f(t,p) 


p(to )  =  Po 


if  we  know  the  value  of  pn  =  p(tn),  this  method  consists  in  calculating  the  new  variable  pn+ 1  =  p(tn  +  At) 
with  the  following  procedure  : 

Pn+l  =  Pn  +  At  K  (3.15) 

K  is  an  evaluated  value  of  the  slope  of  p(t).  For  4th  order  Runge-Kutta,  we  calculate  K  with  a  four-step 
process  : 

K  =  —  (fci  +  2&2  +  2&3  +  Aq)  (3.16) 
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This  formula  introduces  four  coefficients,  Aq,  k2,  and  k 4  defined  as  following  : 


h 

—  f(tm  Pn) 

(3.17) 

k2 

=  f(tn+^Pn  +  -yfci) 

(3.18) 

k3 

At 

=  f{tn+\iPn  + 

(3.19) 

&4 

=  / (^n+1 5  Pn  T  ^3) 

(3.20) 

k\  is  the  slope  of  p(t) 
k2  is  the  slope  of  p(t) 

ks  is  the  slope  of  p(t) 

Pn+\' 

&4  is  the  slope  of  p(t)  evaluated  at  the  end  of  the  interval  (£n+i),  using  the  value  of  k%  to  calculate  pn+ 


evaluated  at  the  beginning  of  the  interval  (tn). 

evaluated  at  the  middle  of  the  interval  (tn+ 1),  using  the  value  of  k\  to  calculate 
evaluated  at  the  middle  of  the  interval  (tn+ 1),  using  the  value  of  k 2  to  calculate 


3.2.2  Spatial  discretization 

To  solve  the  equations  in  space  there  are  three  different  schemes  available.  Two  of  them  are  centered  (4th 
and  6th  order)  and  the  third  is  a  6th  order  compact  scheme.  As  an  example  we  assume  that  (p  =  tp(x), 
and  we  are  estimating  the  value  of  This  description  is  relevant  for  computational  domains  with 
uniform  grids. 

For  the  spatial  derivatives  using  centered  scheme,  the  grid  spacing  Ax  has  to  be  constant.  The 
discretized  equations  are: 


4th  order 


6th  order 


dp  (p{xi-2)  ~  8ip(xi~  1)  +  8ip(xi+i)  -  p{xi+2) 

dx[Xi,~  12  Ax 


(3.21) 


dip  3)  +  9 <p(xi-2)  -  45</?(xj_i)  +  4Sp(xi+l)  -  9 <p(xi+2)  +  <p(xi+3) 

dx{Xi)~  60  Ax 


(3.22) 


When  the  values  of  (p  for  adjacent  points  are  not  available,  such  as  in  the  limits  of  the  domain  for  non¬ 
periodic  boundary  conditions  the  spatial  derivatives  need  to  be  computed  with  a  uncentered  formulation. 

For  the  compact  scheme  the  following  method  described  by  Lele  [5]  is  used: 


ai<Pi~  1  +  bip'i  +  CiPi+i  —  di  (3.23) 

where  oq,  Ci  are  constant  coefficients  and  di  is  a  variable  that  depends  on  the  value  of  ip  at  point  xi 
and  its  neighbors.  Several  compact  schemes  can  be  defined  through  the  choice  of  the  coefficients  oq, 

Q,  and  the  variable  di  [5].  To  achieve  a  6th  order  scheme,  however,  the  following  values  are  used: 


6th  order 
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Similar  to  the  centered  scheme,  the  compact  scheme  uses  values  from  adjacent  nodes.  When  these 
nodes  are  not  available,  for  instance  at  the  boundaries  of  the  domain,  an  alternative  formulation  is 
specified. 
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Chapter  4 

Results  and  discussion 

4.1  Comparison  of  the  thermodynamic  modules 


As  stated  previously,  the  code  YWC  gives  the  user  the  possibility  to  choose  from  various  thermodynamic 
modules  depending  on  the  requirements  of  the  problem.  The  modules  to  be  considered  for  the  treatment 
of  problems  involving  real  gas  thermodynamics  are  REGATH  and  TRC.  It  is  important  to  know  the 
capability  of  each  of  these  modules  to  represent  as  accurately  as  possible  the  properties  of  a  fluid  at 
transcritical  conditions.  Therefore,  in  order  to  compare  the  different  modules,  a  plot  of  density  as  a 
function  of  temperature  for  nitrogen  at  a  given  pressure  was  used.  The  reference  data  was  obtained  from 
the  online  NIST  Standard  Reference  Database  [6]. 


Temperature  [K] 

Figure  4.1:  Nitrogen  thermodynamic  values  obtained  for  a  pressure  of  3.56  MPa  (Pr  —  1.05)  with  the 
two  real  gas  thermodynamic  modules  of  the  code  YWC  and  the  values  from  the  NIST  Standard  Reference 
Database. 


It  is  clear  in  Fig.  4.1  that  the  module  REGATH  starts  departing  from  the  actual  transcritical  behavior 
of  nitrogen  when  the  temperature  decreases  below  200  K.  The  module  TRC,  created  specifically  for  this 
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kind  of  high  pressure,  transcritical  environments,  follows  the  reference  curve  until  120  K,  where  it  starts 
to  deviate.  This  plot  is  useful  to  know  which  ranges  of  temperatures  are  the  most  accurate  for  these  two 
modules  when  the  pressure  of  the  system  is  near  the  critical  point. 


4.2  Assessment  of  the  order  of  the  numerical  scheme 


In  order  to  test  the  accuracy  of  the  numerical  scheme  used,  a  series  of  one-dimensional  simulations  were 
performed.  The  tests  consisted  of  convecting  a  gaussian-shaped  temperature  variation  in  space  through  a 
distance  equivalent  to  100  domain  lengths.  To  achieve  this,  a  constant  velocity  was  imposed  and  periodic 
boundary  conditions  were  used.  In  principle,  if  there  are  no  dissipative  mechanisms  (for  instance,  thermal 
conductivity),  the  two  profiles  (the  initial  and  the  convected  one)  should  always  be  the  same.  However, 
numerical  errors  that  are  proportional  to  the  order  of  the  computational  scheme  appear  in  the  calculations 
and  they  cause  the  profiles  to  differ.  If  the  numerical  error  is  the  only  dissipation  present,  an  improvement 
in  the  error  should  be  achieved  by  increasing  the  number  of  computational  points  in  the  domain,  which 
is  the  equivalent  to  reducing  the  cell  size.  Thus,  the  shape  of  the  gaussian  pulse  after  100  “turns”  was 
compared  to  the  initial  gaussian  shape  and  the  L2  error  was  computed.  For  the  error  the  following 
expression  was  used: 


L2  Error 


\ 


1  Nx 

w  -  w 

X  1=1 


(4.1) 


where  Nx  is  the  number  of  points  that  were  used  in  the  simulation  and  Toi  and  T ^  are  the  tempera¬ 
tures  at  point  i  at  the  initial  condition  and  after  100  “turns”,  respectively.  Fig.  4.2  shows  one  of  the  cases 
with  the  temperature  profiles  at  the  beginning  and  at  the  end  of  the  simulation.  The  initial  condition  is 
presented  in  diamonds  and  the  crosses  represent  the  profile  after  100  “turns”.  The  small  cell  size  (0.0025 
m)  helps  achieve  a  more  accurate  response.  Some  slight  differences  between  profiles  can  be  observed  in 
the  128  to  132  K  range. 


Distance  [m] 


Figure  4.2:  Plot  reflecting  the  evolution  of  a  Gaussian  temperature  profile  of  nitrogen  fluid  after  being 
convected  100  domain  lengths.  The  full  domain  of  the  pulse  (not  shown  here)  is  0.8  meters 
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After  several  similar  one-dimensional  simulations  were  performed  the  L2  error  as  a  function  of  cell  size 
was  retrieved  for  different  scenarios.  In  a  first  assessment,  the  module  AIR  was  used  with  and  without 
thermal  conductivity.  One  can  observe  in  Fig.  4. 3. a  that,  after  a  certain  threshold,  having  a  transport 
property  active  such  as  thermal  conductivity  will  generate  a  numerical  error  that  does  not  improve  by 
decreasing  cell  size.  On  the  other  hand,  when  no  thermal  conductivity  is  present  the  value  of  the  slope  of 
the  error  as  a  function  of  cell  size  is  closer  to  the  order  of  the  numerical  scheme  used  in  the  simulations. 
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Figure  4.3:  Plot  of  the  numerical  error  as  a  function  of  cell  size  in  logarithmic  scale.  The  slope  of  the 
curve  should  be  close  to  the  order  of  the  numerical  scheme  used  in  the  simulation,  (a)  Error  using  the 
AIR  thermodynamic  module  of  YWC.  (b)  Error  using  the  TRC  thermodynamic  module  of  YWC  with 
nitrogen  as  fluid. 


In  a  second  study,  the  same  procedure  was  followed.  This  time  it  was  the  module  TRC  that  was 
analyzed  to  look  for  behavior  resembling  the  order  of  the  numerical  scheme.  The  first  data  points  did 
not  give  results  close  to  the  order  of  the  scheme  and  appeared  as  if  belonging  to  a  fourth  order  scheme. 
However,  after  carefully  reviewing  points  with  smaller  cell  sizes,  a  trend  pointing  to  a  sixth  order  scheme 
was  found  as  seen  in  Fig.  4.3.b. 


4.3  Two-dimensional  direct  numerical  simulations 


The  last  problem  considered  in  this  study  of  transcritical  fluids  using  direct  numerical  simulation  concerns 
two-dimensional  simulations  of  flows  travelling  in  opposite  directions  and  forming  a  shear  layer  in  a 
rectangular  domain.  A  hyperbolic  tangent  profile  is  used  to  define  a  smooth  transition  between  the  flows 
and  avoid  numerical  artifacts  at  the  interface.  A  small  perturbation  is  introduced  in  the  flow  to  initiate 
the  development  of  vortical  structures.  A  diagram  showing  the  main  flow  features  of  this  simulation  is 
shown  in  Fig.  4.4 

The  thermodynamic  module  chosen  for  this  problem  was  REGATH  since  the  conditions  are  very  close 
to  ambient.  The  fluid  was  initialized  at  atmospheric  pressure  with  a  relative  velocity  of  20  m/s  and  a  10 
K  difference  between  them.  The  top  and  bottom  boundaries  were  specified  as  wall  boundary  conditions 
while  the  side  boundaries  were  periodic  boundaries.  The  flow  features  a  perturbation  quantity  in  the 
vertical  direction  with  a  maximum  amplitude  of  2.5  m/s.  The  physical  domain  is  0.2  m  by  0.2  m  and 
the  cell  size  is  1  mm.  The  CFL  number  is  0.7.  Solutions  showing  the  density  field  at  different  times  are 
presented  in  Fig.  4.5. 
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Figure  4.4:  Problem  modeled  using  the  YWC  direct  numerical  simulation  code. 


One  observes  in  Fig.  4. 5. a  the  initial  condition  of  the  flow,  with  two  different  uniform  densities 
throughout  most  of  the  domain.  At  the  interface  there  is  a  layer  which  constitutes  the  transition  between 
both  flows;  this  layer  is  perturbed  by  an  imposed  sinusoidal  vertical  velocity  component  that  creates 
the  characteristic  wave-like  appearance  in  Fig  4.5.b  which  later  evolves  into  a  wave  pattern  with  more 
pronounced  crests  in  Fig  4.5. c.  This  mixing  mechanism  is  then  set  in  motion  resulting  in  the  development 
of  vortical  structures  as  can  be  observed  in  Fig.  4.5.d  which  then  eventually  undergo  vortex  pairing  and 
become  two  larger  vortices  between  50  and  70  ms  as  shown  Fig.  4.6  and  finally  a  very  large  vortex  starts 
developing  just  after  0.1  s  into  the  simulation  (Fig.  4.7). 

In  the  previous  process,  after  the  vortical  structures  are  formed  in  the  first  step  they  continue  to 
entrain  flow  into  the  mixing  layer  which  grows  in  size  and  becomes  very  thick.  The  four  vortices  reach 
a  stable  condition  for  some  time  until  a  perturbation  starts  growing  in  size  (Fig.  4. 6. a)  and  each  pair 
of  vortices  starts  rotating  into  each  other  and  merging  (Figs.  4.6.b  and  4.6.c).  The  process  ends  when 
two  larger  vortices  are  created  as  seen  in  Fig.  4.6.d.  These  same  four  steps  can  be  seen  for  the  same 
simulation  in  Fig.  4.7  when  the  two  vortices  in  Fig.  4.6.d  start  merging  to  create  a  single  very  large 
vortex. 

A  second  simulation  is  shown  in  Fig.  4.8.  It  has  the  same  setup  (CFL  number,  type  of  boundaries, 
domain  size  and  number  of  nodes)  and  presents  the  same  flow  features  shown  in  Fig.  4.4.  The  only 
exception  being  the  temperature  difference  between  the  two  nitrogen  streams.  In  the  first  simulation  the 
difference  is  10  K  and  in  this  second  simulation  the  difference  is  200K.  The  density  ratio  between  the 
denser  flow  and  the  lighter  flow  is  now  2.  Fig.  4.8  shows  density  contours  of  this  simulation  at  the  same 
physical  times  as  the  density  contours  shown  in  Fig.  4.5.  The  processes  taking  place  at  those  times  in  the 
shear  layer  are  practically  the  same  even  though  the  density  ratio  in  the  first  case  is  only  1.03.  Further 
simulations  at  higher  pressures  and  larger  density  ratios  between  the  fluid  streams  should  provide  more 
information  into  the  mechanics  of  shear  layers  under  transcritical  conditions. 
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Figure  4.5:  DNS  results.  Density  contours  at  different  times  for  a  nitrogen  mixing  layer  problem  at 
ambient  conditions,  (a)  0  s,  (b)  0.954ms,  (c)  2.04  ms,  (d)  3.25  ms.  The  plots  correspond  to  a  0.20  m  by 
0.20  m  spatial  domain. 
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Figure  4.6:  DNS  results.  Density  contours  at  different  times  for  a  nitrogen  mixing  layer  problem  at 
ambient  conditions,  (a)  57.0  ms,  (b)  60.4ms,  (c)  63.5  ms,  (d)  66.7  ms.  The  plots  correspond  to  a  0.20  m 
by  0.20  m  spatial  domain. 
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Figure  4.7:  DNS  results.  Density  contours  at  different  times  for  a  nitrogen  mixing  layer  problem  at 
ambient  conditions,  (a)  94.9  ms,  (b)  112  ms,  (c)  122  ms,  (d)  130  ms.  The  plots  correspond  to  a  0.20  m 
by  0.20  m  spatial  domain. 
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Figure  4.8:  DNS  results.  Density  contours  at  different  times  for  a  nitrogen  mixing  layer  with  a  temperature 
difference  of  200K  between  the  two  nitrogen  streams,  (a)  0  s,  (b)  0.954ms,  (c)  2.02  ms,  (d)  3.31  ms.  The 
plots  correspond  to  a  0.20  m  by  0.20  m  spatial  domain. 
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Chapter  5 

Conclusions 


Different  aspects  of  the  direct  numerical  simulation  of  fluid  problems  were  studied  in  this  part  of  the 
report.  This  has  been  accomplished  with  the  help  of  a  new  modeling  tool  developed  at  the  EM2C 
laboratory.  The  code,  YMC,  has  been  developed  with  a  large  variety  of  applications  in  mind  and  the 
need  to  understand  many  of  its  functionalities  was  emphasized.  In  particular,  the  use  of  thermodynamic 
modules  to  address  different  types  of  problems  led  to  the  analysis  of  the  equation  of  state  of  the  two 
modules  that  dealt  with  real  gas  formulations,  which  is  the  type  of  thermodynamic  capability  needed 
for  the  larger  goal  of  the  current  research  efforts.  Also,  the  results  of  one-dimensional  validation  cases 
were  used  to  confirm  the  accuracy  of  the  numerical  scheme  to  sixth  order,  which  was  the  order  of 
the  scheme  employed  to  solve  the  governing  equations.  Furthermore,  two-dimensional  simulations  of  a 
nitrogen  mixing  layer  were  performed  to  analyze  the  behavior  of  the  flow  in  the  interface  and  study  its 
behavior.  The  different  tasks  performed  during  this  research  helped  understand  and  assess  the  different 
characteristics  of  direct  numerical  simulation  methods,  and  later  apply  the  obtained  knowledge  to  the 
solution  of  a  relevant  problem. 
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N  omenclat  ur  e 


CFL 

CNRS 

^ij 

Dk 

DNS 

e 

EM2C 

FSA 

LES 

P 

NIST 

P,P 

Pc 

Pr 

Qi 

P 

Re 

REGATH 

SGS 

T 

t 

T 

TRC 

Ui 

Uk 

Xi 

Ax 

Yk 

YWC 


Courant-Friedrichs-Lewy  (a  condition  of  convergence  in  partial  differential  equations) 
Centre  National  de  la  Recherche  Scientifique  (France) 

Kronecker  delta  (function  that  can  be  only  0  or  1  depending  on  the  values  of  i  and  j ) 
diffusivity  of  species  k 
Direct  Numerical  Simulation 

specific  energy  or  total  energy  divided  by  mass  (e  =  E/m) 

Energetique,  Moleculaire,  Macroscopique  et  Combustion 
Faculte  des  Sciences  Appliquees 
Large  Eddy  Simulation 
dynamic  viscosity 

National  Institute  of  Standards  and  Technology 

pressure 

critical  pressure 

reduced  pressure  (Pr  =  P/Pc) 

heat  flux  vector 

density 

Reynolds  number  (Re  =  pul / p) 

REal  GAs  THer  mo  dynamics 
Subgrid  Scale 
temperature 
time 

viscosity  stress  tensor 
TRansCritical 

abbreviated  notation  for  the  three  velocity  components  u,  v  and  w 
net  production  rate  of  species  k 

abbreviated  notation  for  the  three  Cartesian  coordinates  x,  y  and  z 
length  between  nodes  in  a  computer  simulation,  cell  size 
mass  fraction  of  species  k 
Yes  We  Can  (an  EM2C  DNS  code) 


101 


